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Abstract 

Generalized Estimation Equations (GEE) are a well-known method for the analy¬ 
sis of categorical longitudinal responses. GEE method has computational simplicity 
and population parameter interpretation. In the presence of missing data it is only 
valid under the strong assumption of missing completely at random. A doubly ro¬ 
bust estimator (DRGEE) for correlated ordinal longitudinal data is a nice approach 
for handling intermittently missing response and covariate under the MAR mecha¬ 
nism. Independent working correlation is the standard way in DRGEE. However, 
when covariate is not time stationary, efficiency can be gained using a structured 
association. The goal of this paper is to extend the DRGEE estimator to allow 
modeling the association structure by means of either the correlation coefficient or 
local odds ratio. Simulation results revealed better performance of the local odds 
ratio parametrization, specially for small samples. The method is applied to a data 
set related to Rheumatic Mitral Stenosis. 

Keywords: Correlation coefficient; Doubly robust estimators; Generalized estimating 
equations; Local odds ratio; Missing at random. 


1 Introduction 


Longitudinal data arise when each individual is measured repeatedly through time. 
These repeated responses form a cluster and it is expected the response within each cluster 
to be correlated. A popular approach for the analysis of longitudinal data is the Gener¬ 
alized Estimating Equation (GEE) method, proposed bv lLiang &: Zegerl (119861) . The goal 
of this proce dure is to estirnate hx ed parameters without specifying a joint distribution 
for the data ( Nooraee et all 2014 ). GEE method requires only the correct specihcation 
of the response’s mean structure for the parameter estimator to be consistent and asymp¬ 
totically normal. The attractive feature of GEE is that the association parameters among 
repeated measures are taken as a ‘nuisance’ parameters and, unlike maximum likelihood 
methods, mean parameter estimates are not sensible to the specihcation of the association 
structure. Furthermore, GEE method allows marginal interpretation of the parameter of 
interest and it has computational simplicity. 
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A simple way of analyzing such correlated data is to consider an independence working 
assumption for the repeated responses. However, when the covariat e design is not time - 
stationary it will lead to inefficient marginal regression estimates f Lipsitz et all 1994 1. 
In the presence of time-varying covariates efficiency can be gained assuming a different 
correlation structure. Nevertheless, modeling the association in ordinal data is not a 
simple task. Different approach es have been propo sed to estimate the association pa- 

) provided moment estimators for a 
teOOGll proposed an approach which 


rameters for ordinal responses. iLipsitz et an 11199 


variety of correlation matrices, while iParsons et al 


estimates the correlation vector by minimizing the logarithm of the det erminant of the 
covariance matrix of the hxed parameters. Instead of using correlations, Lnmlevl ( 1996 1 
proposed using a c ommo n global odds to reduce the number of a ssociation paramet ers. 


Heagertv fc Zegerl (119961) extended the ALR method proposed by ICarev et al\ (119931) to 


ordinal re sponses using a s econd set of estimating equations for the global odds ratio. 
Recently, iTouloumis et al\ (1201311 considered a family of association models to estimate 
local odds ratios as a measure of ass ociation. A comparis on study of different working 
association structures can be found in Nooraee et all ( 2014 1. 

In the presence of missing data, inferences are valid if the missing ness mechanism 
is missing completely at random (MCAR), as dehned bv iRnbinI (1197611. When d ata is 
MAR, one can adopt mult iple imputation GE E (MIGEE) (jhittle fc Rubin . 1987 1 or a 


weighted version (WGEE) ( Robins et al . 1995ll . These single robust versions of GEE for 


incomplete data require the correct specihcation of the weight model of GEE fWGEEl or 


the imputation mode l fMIGEEl. Doub l y robust estimators tDRGE El (IGarpenter et al. 


( 2006 1. TsiatisI ( 20061) . Seaman fc Gopas ( 2009ll . Ghen fc ZhoiJ ( 201lll ) combine ideas from 
these two approaches. For consistency, it requires only the weight or the imputation model 
to be correctly specihed, providing more flexibility for the modeler. 

This work was motivated by the Rheumatic Mitral Stenosis study in which a cohort 
of 164 patients with rheumatic mitral stenosis (a narrowing of the mitral valve in the 
heart) were referred for treatment at Hospital das Glinicas of the Federal University of 
Minas Gerais, Brazil. The response of interest was the functional classihcation (NYHA), 
major determinant of quality of life and survival of the individual. The main objective of 
that study was to characterize the improvement of the functional classihcation over time. 
This study was characterized by an arbitrary pattern of missing data. Response and 
a particular covariate (atrial compliance) were missing for some patients and the MAR 
mechanism seems to be a reasonable assumption for this data. 

In the current paper, we consider a doubly robust approach for the analysis of longitu¬ 
dinal ordinal data with intermittently missing response along with a key covariate that is 
MAR. The cumulative logit model was used for the marginal means. We extend the dou¬ 
bly robust estimator to ac commodate two par ametrizations of the assoc iation structures : 
the c orrelation coefficient ( Lipsitz et ajl 1994f ) and the local odds ratio ( Tonlonmis et al . 
2 OI 3 I) . Efficiency and accuracy of the proposed estimator are compared under these asso¬ 


ciation structures. 

The paper is organized as follows. In Section [2] are dehned the notation for GEE with 
fully observed data and discuss the two parametrizations of the association structure. Sec¬ 
tion |3] outlines WGEE and MIGEE approaches. The proposed methodology is established 
in Section 01 A simulation study is presented in Section |5l in which the hnite-sample 
biases and mainly standard errors are compared for the standard GEE, MIGEE, WGEE 
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and doubly robust versions, under both the correlation and local odds parametrizations. 
Data arising from the Rheumatic Mitral Stenosis study are analyzed in Section [6l Paper 
ends with a discussion and future directions in Section [71 


2 Notation and GEE for Complete Data 

In this section it is introduced the generalized estimating equations for the analysis of 
fully observed ordinal data. Section 12.11 establishes the model and notation for longitu¬ 
dinal ordinal data. Section 12.21 presents two competing ways of modeling the association 
structure in GEE. 


2.1 GEE for Longitudinal Ordinal Response 

Let Oit G {1, 2,..., J} be the ordinal response for i-th. subject (i = 1,..., n) at time 
f (f = 1,..., Tj, Ti <T). As the response has J levels it can be dehned as Y^j = /{On = j) 
for j = 1,..., J, where I {A) denotes the indicator function. Yuj is converted into the 
equivalent (J — l)-variate vector Ya = {Ym ,..., and let = (Ejf,..., Y^J'^ 

be the stacked response vector. When J = 2 the response is binary and Y^is a. scalar. Let 
Xi denotes the time-stationary covariate for the i-th subject, and Zi = ..., Zjj,,)'^ 

a Ti X q matrix of explanatory variables. 

The marginal distribution of Yn is assumed to be multinomial = 1); 

is 

J 

i=i 

where finj = = E{Yitj\Xn Z^ (3) = Pr{Oit = j\XnZi,f3), is the probability of 

response j at time t and /B is a p x 1 vector of parameters. Two common choices for 
modeling pnj are the cumulative logit and probit models. In this work it is assumed a 
cumulative logit link, that is. 


logit [Pr{Oit < j\Xn Zu)] = /Soj + Xi(]^ -F Z^/3^, j = 1,..., J - 1. 


( 2 ) 


Eormulation in ([2]) implies a proportional odds model flMcCullaghl . 119801 ). In such 
model the interpretation of /3 is the same regardless of the number of categories (i.e., it is 
invariant to combination of categories). A desired featur e is that the exponential of the 
parameters is interpreted as an odds ratio f Agrestil . 2013). 

The main interest is to make inferences related to the regression parameters 
/3 = (/5oi, ■ ■ ■ 5 /9o,j-i, 0^)'^ associated to the (J — 1) x 1 marginal probability vectors 

E{Yit\Xi^ Zi) {pitii • • • ■ 

Hn is grouped to form a vector E{Yi\Xi, Zi) = pLi = , pijrpj'^ with the same 

dimension of 

_ In order to estimate 0 , generalized estimation equations are used (ILiang fc Zeger 


(1198611 : iLipsitz et al\ (119941 ). iToulonmis et al\ (120131 )). which takes the form 


5^C/.(/3,a) = 5^DW-‘(V-M.) = 0, (3) 

i=l i=l 
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where Di = and Vi = Vi{^, a) is a Tj( J—1) xTi( J—1) covariance matrix for Yi. The 
vector parameter a expresses a ‘working’ assumption about the correlation/association 
structure. 

Under mild regularity conditions, correct specihcation of the marginal mean model in 
(j2]), and provided that a y^-consistent of a is available, Liang and Zeger (1986) proved 
that the estimator obtained by solving ([S]), is consistent and \/n(/3 — /3) converges 
in distribution to a multivariate normal distribution with mean vector 0 and covariance 
matrix 

Vp = lim nSo'SiSo \ (4) 

n—>-oo 


where Sq = Eti D^Vr^DJ, and = j:i^D,V-^CoY{Y,)Vr^Df. 

In practice, the “sandwich” covariance matrix in (jl]) is calculated by ignoring the 
l imit and replacing (13, 0 .) and Cov(l^j) by and {Yi — p,j)(Yi — p^i)'^, respectively 


f Toulonmis et ai . 20131 ). 


2.2 Estimation of the nuisance parameter vector and covariance 
matrix 

The ‘working’ assumption is ‘independence’ when no correlation is assumed between 
pairs of the response of each individual, that is, 0 = 0. Independent Estimating Equations 
(lEE) are proved to be efficient only when covariates are constant over time or if the 
independence structure is actually true f Lipsitz et all 1994 1. In this case, score equations 
of the ML method for the regression vector (3 are identical to the lEE if all observations are 
treated as independent. On the other hand, when there exist within-individual association 
or time-varying covariates, efficiency can be gained by modeling the correlation structure. 
Hence, a number of proposals have been formulated to estimate a. These alternatives 
dif fer in the effic i ency o f estimating the covariance matrix and computational simplicity. 

Lipsitz et al. f 1994fl dehned a as the correlation coefficient, and sug gested the use 


of the method of moments to estimate a number correlation structures. Parsons et al. 


f 2nn6fl proposed an approach t hat estimate th e correlation vector a by minimizing an 


objective function Q(a \0, y). iLumleyl (119961 ) proposed using a common global od ds. 
Heagertv &: Zegerl (119961 1 extended the ALR method proposed by ICarev et all (1199311 to 


ordinal responses by using a se cond set of estimating equations for the global odds ratio. 
Finally. iTonlonmis et al\ (120131 ) identihes ck as a ‘nuisance’ parameter vector that contains 
the marginalized local odds ratios structure. They employed a family of association models 
in order to develop meaningful structures for the ordinal response. 


2.2.1 Correlation Coefficient 


Lipsitz et al. ( 1994{ l suggested a method that constrains the correlations at different 
times between two categories of the response. In their approach the weight matrix Vi 
is decomposed into the form Vi{(3,ct) = Fy^{(3)Ci{ct)Fy^{f3), where Fi is a matrix 
containing marginal variances, F^, given by 

Fit diag [/ijfi(l /^iil); ■ ■ ■ 1 —lip l)] ; 

and Ci is equal to the marginal correlation matrix. The unknown elements of Ciicx) are 
the elements of the (J — 1) x (J — 1) matrix p^^^,icx). Thus, is parametrized as 
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Pitt' = Corr{Yu,Yit'). The (J - 1) x (J - 1) diagonal blocks of Vi are F , 
with Vit = diag^Hit) — (J — 1) x (J — 1) off-diagonal blocks of Ci{a.) are 

Pitt' = Pitt'{^)i which represents the correlation between Yu and Yui, t ^ t'. The vector 
a is a parameter vector associated with the model for 
Dehne the Pearson residual vector, eu as 

ea = F-P' (Y„ - . 


Then it follows that 


Citt'{a) = Corr{Yu,Yit>) = E^eueJu)- 


In order to reduce the dimension of the correlation vector, Ihipsitz et all (119941) assumed 
an uniform correlation structure over the individuals. The use of the method of moments 
was suggested for a variety of correlation matrices such as 

• exchangeable: Pnu = p, for all t < t'; 

• 1 -dependent: Put+i — Pt) ^ = 1) • • • > ^ “ 1) cind pm' = 0 otherwise; 

• handed: p^i = Pr, when |T — t| = r, for r = 1,..., T — 1; 

• unstructured: put/ = Put'- 

The estimate a is plugged into ([S]) and a solution is found for /3. The solution might 
be obtained by a Fisher-scoring algorithm. 

The correlation coefficient parametrization ignores the scale of the response variable 
and may resul t in loss of i nformation regardin g the correlation between the variables 
f Lumlev . 1996h . Moreover, Lipsitz et al. f 1994tl observed that the ‘working’ correlation 
matrix is not always positive definite, which may result in a breakdown of the Fisher 
scoring method. This is specially true over unstructured correlations matrices and small 
sample sizes. When the given model for Puti contains too many parameters, the resulting 
estimates of (3 and a may be highly variable. 


2.2.2 Local Odds Ratio 


Inst ead of modeling the corr elation coefficients, some authors (see, for example. iLnmley 


(119961) . iTouloumis et al\ (120131) ) modeled the off-diagonal Vat' through the joint probabil¬ 
ity of the responses Yu and Yu/, for t ^ t'. The covariance of Yuj and Yuj/ can be writen 


as 


CoviXitj, Yu'f\Xi, Zu, (3) 


EiXuj, Yu'j/\Xi, Zit, (3) — E{Yuj\Xi, Zu, f3)E{Yitiji\Xi, Zu, (3) 

f^itj f^it^ ' 


The product of the first moments can be calculated through the specification of the 
marginal model ([2]). The joint probabilities pujt'j' = P(Xitj = '^,Yit'j'\Xi, Zu, (3) can 
be modeled by an association vector that descr ibes the associat i on str ucture for Vi = 
1,..., n, t 7 ^ t' = 1,..., Tj and j,j' = 1,..., J ( Tonlonmis et a,l. . 20131 ). Estimation of 
the se joint probab ilitie s can be achieved t hroug h the global odds rati o (see, for exam - 


ple, iLumlevI (119961) and iHeagertv &: Zegerl (ll996f) ) or a local odds ratio (iTouloumis et al 


5 







































20131). iTouloiimis et al\ (120131 ) argue that the local odds ratio provides the best parametriza- 


tion because it is variation independent to (3 and they produce valid and unique positive 
joint probabilities. 

Let’s introduce the local odds ratio of Touloumis et al\ f 2013h . Denote the L = T(T — 
l)/2 time pairs as (1, 2), (1, 3),..., (T — 1, T), where T = max{Ti ,..., T„}. Let Ftjt'j' = 
j' the observed frequency of the cutpoint (j, j') at the (t,t') time pair of the 
marginalized table. Dehne as the local odds, that is, 




Ftjt'j',j'+i 

Ftj+i,t'j'Ftj^t',j'+i 


for j,j' = 1,... , J — 1, and let a be the L x {J — 1)^ vector consisting of the local odds 
ratio 

Ot = (6*1121, • • • , 6*ii2(J_1), . . . , d(^T-l)lTl, ■ ■ ■ , 9{T-1){J-1)T{J-1))'^■ 

The association vector a can be estimated by htting a loglinear model for the counts 
{Ftjt'j'} simultaneously to all possi ble L marginalized cont ingency tables and then calcu¬ 
lating the implied local odds ratio (ITonlonmis et oi.l . 120131 1. For notational reasons, let A 
and B be the row and and column variable, respectively, and let G be the group variable 
with levels being the L ordered pairs. Assuming a P oisson sampling scheme to the L sets 
of J X J contingency tables, £t the RC type model flBecker fc Cloggl . Il9891 


'I - \ I I I \G I \AG I ^BG I (t,t') (t,t') (t,t') 

^Og[jtjt'j'> — A + Aj -h Aj, + Af^tt') + A(tt') + A(^tt') + ^j ^j> ) 


(5) 


where ^ • J = 1) • • • > <^1 Eire the score parameters for the J response categories at 

time pair (f, F), and {ftjt'j' '■ j,j' = ■ ■ ■ y J} are the expected frequencies. The maximum 

likelihood estimate of a are obtained by treating the L marginalized contingency tables as 
independent. By imposing identihability constraints on the regression parameters in (]5|), 
the resulting unrestricted local odds ratio are determined by \og{6tj't'j') = ^ — 

^ ~^j'+i)y where the intrinsic parameter (pd.*') measures the average association 
of the marginalized contingency table. To increase parsimony, common unit-spaced score 
parameters {Uj ’ = j) are usually assumed. Main options for the marginalized local odds 

ratio structures include 

• uniform: \og{9tj>t'j>) = estimates a single parameter; 

• time exchangeability: \og{9tj't'j') = estimates L parameters; 

• category exchangeability: \og{9tj't'j') = g^iyj — r'j+i){r'j> — i^y+i); estimates J — 1 
parameters, and 

• unstructured: \og{9tj't'j') = (p9A)(^i/9 ’^) _ — that requires L(J —1) 

parameters. 

Conditional on ci, and the marginal specification ([2]), the joint probabilities fiujt'j' 
are estimated based on the adopted local odds ratio structur e using the IPFP (Iterativ e 
Proportional Fitting Procedure). This algorithm, proposed by ( Deming &: Stephanl . 1940 ). 
is used to obtain jiujt'j' through the marginals fittj and fiu'j'- iTouloumis et al\ fl2013l) 
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proved that the IPFP solution preserves local odds ratios of the initial values as long as 
they are positive. Hence, it is straightforward to calculate the weight matrix Vi and the 
estimating equations in ([3]) can be solved with respect to (3. 

An advantage of the local odds ratio parametrization over the correlation is that the 
local odds ratio and the marginal regression vector are variation independent. This means 
that /3 estimates are less sensible to a possibly wrong specihcation of a. As opposed to 
the correlation parametrization, the estimation of the association parameters does not 
depend on covariates and, as long as it is based on maximum likelihood models, the cx. 
estimates are consistent under MAR. Thus, no adjustment need to be done on o: obtained 
with available data. 


3 Available Approaches for Missing Data 

Section 13.11 presents a series of assumptions related to mechanism causing data to be 
missing and necessary to be considered in order to build valid estimators. Multiple impu¬ 
tation and weighed generalized estimation equations are two commonly methods available 
for missing data under MAR mechanism. These methods are presented in Sections 13.21 
and 13.3[ respectively. They serve as the basis for the construction of the doubly robust 
estimator, presented in Section HI 

3.1 Missing Data Framework 

In this work it will be assumed that the time-stationary covariate X* may be missing 
for some subjects whereas the explanatory variables Zi are fully observed. 

For each occasion t it can be dehned Ru = 0 if On and Xj are missing, Rn = 1 if On 
is missing and X* is observed, Rn = 2 if On is observed and Xj is missing, and Rn = 3 if 
On and Xj are both observed. Let Ri = {Rn ,..., RiTi)'^, and Rn = {Rn, ■ ■ ■, Ri^t-i)- 
The marginal probability Pr{Ri = ri\Oi,Zi) can be obtained through conditional 
models of the form Pr{Rn = rn \Rit, Oi, Zi). This general formulation encompasses 
MCAR, MAR and MNAR mechanisms. In particular, the MAR mechanism requires 

Pr{Ri = ri\Oi, Zi) = Pr{Ri, = ri\0°, Zi), (6) 

where 0° denotes the observed components of Oj. The following natural further assump¬ 
tion is made 

Pr{Rit = rit\Rit, Oi, Zi) = Pr{Rn = rn\Rit, On, Z°^), (7) 

for each time t, where 0°^ and Z°^ are the histories of observed responses and covariates 
up to time t — 1. 

Let T^n = Pr{Rn = 3|Oj, Zi) be the marginal probability of observing both Oi and Xj 
at time t, given the entire vectors of responses and covariates. Then, Tin is expressed by 

^ ^ Pr{Rn 3, Rj Tii_\, ..., Rn rn\Oi, Zi). 

This marginal probability can be expressed in terms of the conditional probabilities 
Pr{Rn = k\Rn,Oi,Z i), for k = 0,1,2,3. Throughout this paper it will be required 
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the so-called positivity assumption, that is, nu must be bounded away from zero. This 
condition is needed in order to guarantee the existence of -^n-consistent estimators of (3 
( Robins et al . 1995 1. 


3.2 Multiple Imputation Generalized Estimating Equations 

A imputation model commonly used to handl e intermittently rn i ssing r esponse and 
covari ate, is imputation using chained equations ( van Buuren et ail ( 1999 1. van Buuren 
( 2007ll l. which is more commonly referred to as full conditional specihcation (FCS). This 
approach specihes conditional distributions for each incomplete variable, conditional on 
all others variables in the imputation model. Starting from an initial imputation, FCS 
draws imputations by iterating over the conditional densities. 

Denote by f3^ and Um, respectively, the estimate of f3 and its covariance matr ix frorn 
the G EE analysis of the m-th completed data set, (m = 1,... ,M). Following (iRubin . 
198711 , the combined point estimate for the parameter of interest (3 based on MI is simply 
the average of the M complete-data point estimates 

1 ^ 

^MI ~ Pmi 

m=l 


and an estimate of the covariance matrix of i3j^j is given by 

'M + 1' 


Umi = W 


M 


B, 


where 


^ M ^ M 


m=l 


m=l 


3.3 Weighted Generalized Estimating Equations 

Robins et a,l\ ( 19951 ) proposed a class of weighted estimating e quations to allow for 


MAR mechanism. In binary longitudinal data, IChen fc Zhoul (120111 ) extended the method 
to accommodate arbitrary patterns of both missing response and covariate. Their method 
was adapted here for longitudinal ordinal responses. 


Dehne a weight matrix Aj = [5, 




t = 1 , 


= 1 , 


, T,-, where 


Siw = {HRit = 1, Rif = 3) + I{Rit = 3, Rif = 3)} /iTiu' for t 7 ^ t', 6^ = I{Rit = 
and Tiitf = Pr{Rit = l,Rit' = 3\Oi,Zi) + Pr{Rit = 3, Rif = 3\Oi,Zi). In order to 
construct the weight matrix Aj the conditional probabilities Pr{Rit = k\Rit, Oi, Zi) are 
decomposed as the product of two separated logistic models. The hrst one models the 
probability of observing the potentially missing covariate Xi whereas the latter models 
the probability of observing conditional on observed responses and covariates up to 
time t — 1. 

The main idea of weighted generalized estimating equations (WGEE) lays on weighting 
the individual contribution to the estimating equation by introducing the weight matrix 
Aj into the covariance matrix Vi. This task can be accomplished by different ways 
























depending on the parametrization of the association structure. The general WGEE for ^ 
are given by 

n n 

a, V’) = D,M,{Y, - m,) = 0. (8) 

i=\ i=l 

When the correlation coefficient is adopted, the Mi matrix is decomposed as Mi = 
Fi • Aj)E'j where A ■ B = [au ■ bu] denotes the Hadamard product of matrices 

A = [au] and B = \bu]. Under the local odds approach, the matrix Mi is given by 

M, = yri . A,. 

Conditionally on a consistent estimate of the correlation/association structure ck, a 
consistent estimate for (3 can be obtained by solving ([8]), under the correct specihcation 
of the missing data model. 

In the presence of missing data a consistent estimate for the correlation parame¬ 
ters can be obtained by dehning a weighted observed pair of Pearson residual vector 
as Cii = eu{I{Ru = 3)/7iu). Then, a moment-based estimator can be constructed using 
the weighted pair contribution. For instance, an exchangeable correlation estimate can be 
obtained through 


pitt' 




n 

CtCi' • 


t'>t 




It is easy to show that put' is an unbiased estimator. 


4 Doubly Robust GEE for Longitudinal Ordinal Data 


Some authors (e.g., Scharfstein et a,l. f 19991 ). iTsiatisI fenObl )) noted that adding a term 


of expectation zero, say to the inverse probability weighted estimators would still re¬ 
sult in consistent estimates under the MAR mechanism. The solutions of these augmented 
estimating eq uations give ri s e to t he so-called doubly robust estimators. 

Following IChen fc ZhonI (120111) , the optimal <popt for missing response and covariate 


is given by = E(^YY^^x”^iY°,x°,Zi,Ri) {DiNi{Yi - Hi)}, where and X™ denote 
the missing components of Yi and Xj, respectively. When the correlation coefficient 


parametrization is adopted, Ni is dehned as Ni = ^ • (11^ — A*)} F] 

where 1 is a vector of I’s of length Tj(J — 1). With local odds Ni can be dehned as 


n /2 


N, = V 


-1 


lU 


Conditionally on a consistent estimate of the correlation/association structure a, an 
improved estimate for (3 can then be obtained by solving the estimating equations 


n 


i=\ 


[DiMi[Y , - {r>iNi{Y , - Mi)}] — 0- 


(9) 


The estimator for /3 in ([9]) is doubly-robust in the sense that it is consis tent if at least 
one o f the missing data model or the covariate model is correctly specihed (jChen fc Zhoil 

201 ll) . 
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The referred expectation in the second part of (El) is over the conditional distribution 
of (1^™, Xf, Zi, Ri), which can be written as 

p{YT = yT,xr = xT\Yix^,z„R,-f3*,'Y) = p{YT = vT, = xT\Y°, x:^, z,-p* 

= P{YT = yT\YlX, = x„Zf,f3*) 

xP{Xr = xT\YlX^,Z,-j). 

The multivariate distribution P{Y’^ = y'^\Y°,Xi = Xi,Zi;(3*) is expressed by a 
product of univariate ordinal models. When X is discrete, the second term in ([9]) can be 
written as 

E(YY^^XY'\Y°,x°,z„R,){DiNi{Yi - y.^)} = ^ Wi:j;y{DiNi{Yi-y.i)}, 


iyT,^T) 


where the weight Wixy is given by 


W- 


ixy 


= P{YT = yT\YlX, = x,, Zp/3*) X P{X^ = xT\Y^„X°, Z,;j). 


In the case of X continuous, the second term in (El) takes the form 


D,N,{Yi - = f w,xy. {D,Ni{Y, - fi,)} dYTXr, 

This expectation can be cumbersome, depending on the missing data pattern. In such 
case, instead of using numerical integration, a Monte Carlo method can be applied to 
approximate the corresponding integral. 

Inspired by doubly robust ideas, we constructed the following estimator for the corre¬ 
lation structure 


E{Y”^,XY^\Y°,X°,Zi,Ri) { 


n /-, N n 


2=1 t'>t 


2=1 t'>t 






where = Ti{Ti — 1) — p and 0 < a; < 1. 

A sandwich estimator for the standard error of /3 is given in Appendix. 


5 Simulation Study 

In this section a simulation study is presented in order to investigate the performance 
of the proposed method under the two parametrizations of the association vector as well 
as its robustness to misspecihcation of the predictive models. It is considered a study 
with Ti = T = 3 repeated ordinal measures (with three categories) and two covariates 
(one quantitative and other qualitative). The true marginal model is 

logit Pr{Oit < j\Xit, Zu) = f3oj + jdiXi + f32Zit, j = 1, 2. (10) 

where Z^ ~ X(0,1/2) for f = 1, 2, 3. 

The binary covariate X^ may be missing for some subjects and is generated according 
to 

logit Pr{Xi = l\Zii) ='jo + liZn. (11) 
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It is assumed that /Sqi = —0.4, / 5 o 2 = 1-2, /5i = —0.35, 132 = 0.35, 70 = log{l), 
71 =2. The correlated ordinal responses were generated using the NORTA method 


f Anestis Touloumis . 20131) with constant correlation between the latent vectors set equal 
to p = 0.7. 

In order to model Ra two new indicators were dehned. Let the indicator of ob¬ 
serving Xi and Rfj the indicator of observing Ou- The response variable in the hrst time 
occasion was allowed to be fully observed. The model for Rf was dehned as 


\Pr{Rf = 0) 


= ro+rl0^l+r2Zil, 


( 12 ) 


and the model for was taken as 


^^3 ^ = V'o + + ^lI{Rlt_^ = 1) + ^ = 2, 3, (13) 

where 0*^_i = if is observed and 0 otherwise. The true values are taken as 

-00 = 1-2, '01 = —1-5, -02 = —1-5, -00 = 0.6, ip\ = —1.5, -01 = 2.5, and 'ipl = —1.3. It was 
observed about 30% of missing observations under this setup. 

For comparison purposes, it was considered ordinary GEE for the complete and avail¬ 
able data, r espectively, weighted GEE (WGEE), multip le imputation (MIGEE) by chained 
equations ( van Buuren fc Groothuis-Oudshoorn . 201l[l with M = 10 multiple imputa¬ 
tions, and the doubly robust versions (DRGEE). The primary goal of this simulation 
was to compare the performance of the above mentioned methods under the correlation 
and local odds ratio parametrizations. Three correlation structures (independent - ind, 
exchangeable - exch, and unstructured - unst) and four local odds ratio structures (uni¬ 
form - unif, category exchangeability - cat.exch, time exchangeability - time.exch, and 
unstructured - RG) were to be compared. Under independence, the estimates from the 
two parametrization are identical. 

In order to investigate robustness of these methods, the predicted models were also 
misspecihed by omitting the covariate Zi from the covariate model (HI]) or the missing 
data model (|T^ . 

Let S = 1000 the total number of Monte Garlo replications. Whenever, in a given 
iteration, a working association structure (G) failed to converge a new sample data were 


generated. Denote by (3^ the corresponding GEE estimator at the r-th Monte Garlo repli¬ 
cation and let be the arithmetic mean, /S'" = I/S' Pr- To evaluate the consistency 
of the competing methods the relative bias, dehned as 100 x (/S'" — (3)/(3, was calculated 
for each parameter. Interest is in quantifying the gain in efficiency by the working associ¬ 
ation structures over the independence structure. The Monte Garlo relative efficiency was 
dehned as EP{I31)/ 33P{l3^), where EP{l3^) is the standard error of (3^ based 

on the estimated robust covariance matrix under the (G) working association structure. 
The estimated coverage probability for a nominal 95% level based on the asymptotic nor¬ 
mality of the GEE estimators is also reported. Simulation was conducted for sample sizes 
n = 50,150, 300 and 600. Results are presented for sample size n = 300 subjects. 

Table [U presents the simulation results associated with ordinary GEE for available 
data and incorrectly specihed methods, those in which the covariate Zi was omitted 
from their predictive models. For each method the hrst three lines refer to correlation 
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Table 1: Evaluation criteria for misspecified models. Results for n = 300 and S = 1000 
simulations._ 




Relative Bias 


Relative Efficiency 


Empirical Coverage 


Structure 

/3oi 

/3o2 


Z 

Poi 

/3o2 


z 

Poi 

P02 

X 

Z 

ind 

121.1 

-31.5 

-23.9 

60.2 

1.00 

Available 

1.00 1.00 

1.00 

0.15 

0.36 

0.93 

0.73 

exch 

94.0 

-25.2 

-26.4 

32.3 

1.02 

1.01 

1.01 

1.20 

0.35 

0.55 

0.94 

0.86 

unst 

94.3 

-25.9 

-31.4 

26.7 

1.03 

1.01 

1.02 

1.19 

0.31 

0.48 

0.90 

0.90 

unif 

93.7 

-23.9 

-27.8 

26.4 

1.02 

0.99 

1.01 

1.23 

0.33 

0.57 

0.92 

0.88 

cat. exch 

93.8 

-23.9 

-28.0 

25.5 

1.02 

0.99 

1.01 

1.24 

0.33 

0.57 

0.92 

0.89 

time.exch 

93.7 

-24.0 

-27.9 

26.3 

1.02 

0.99 

1.01 

1.24 

0.36 

0.56 

0.93 

0.89 

RC 

93.1 

-24.2 

-30.8 

22.5 

1.02 

1.00 

1.01 

1.24 

0.34 

0.56 

0.92 

0.91 

ind 

19.9 

-5.4 

-28.9 

32.8 

1.00 

WGEE(r-) 
1.00 1.00 

1.00 

0.94 

0.93 

0.92 

0.89 

exch 

17.9 

-5.1 

-29.9 

22.2 

1.00 

0.99 

1.00 

1.17 

0.94 

0.94 

0.94 

0.91 

unst 

21.0 

-6.5 

-35.0 

15.2 

1.01 

0.99 

1.01 

1.16 

0.93 

0.92 

0.91 

0.94 

unif 

18.8 

-5.1 

-28.6 

22.7 

1.00 

0.99 

1.00 

1.16 

0.95 

0.94 

0.92 

0.92 

cat. exch 

19.0 

-5.2 

-28.8 

22.0 

1.00 

0.99 

1.00 

1.17 

0.95 

0.94 

0.92 

0.92 

time.exch 

18.3 

-5.0 

-27.7 

22.8 

1.00 

0.99 

1.00 

1.17 

0.95 

0.93 

0.93 

0.92 

RC 

18.7 

-5.6 

-31.4 

19.8 

1.01 

0.99 

1.00 

1.17 

0.94 

0.93 

0.92 

0.91 

ind 

20.9 

-7.1 

-48.0 

-8.3 

1.00 

MIGEE(x-) 
1.00 1.00 

1.00 

0.93 

0.93 

0.92 

0.94 

exch 

19.9 

-6.8 

-47.9 

-1.6 

0.99 

0.99 

1.00 

1.22 

0.94 

0.94 

0.93 

0.95 

unst 

21.5 

-7.6 

-51.9 

-5.9 

0.99 

1.00 

1.00 

1.22 

0.94 

0.93 

0.92 

0.94 

unif 

19.7 

-6.7 

-46.8 

-2.2 

1.00 

1.00 

1.00 

1.22 

0.95 

0.94 

0.93 

0.94 

cat. exch 

19.9 

-6.8 

-47.0 

-3.1 

1.00 

1.00 

1.00 

1.23 

0.95 

0.94 

0.93 

0.94 

time.exch 

20.0 

-6.8 

-46.6 

-2.8 

1.00 

1.00 

1.00 

1.23 

0.94 

0.94 

0.94 

0.94 

RC 

20.1 

-7.3 

-49.7 

-5.8 

1.00 

1.01 

1.00 

1.23 

0.94 

0.93 

0.92 

0.93 

ind 

17.0 

-5.5 

-41.4 

-8.9 

1.00 

DRGEE( 

1.00 

x ~, r~ 
1.00 

■) 

1.00 

0.93 

0.92 

0.89 

0.94 

exch 

15.7 

-5.1 

-40.7 

-2.4 

0.99 

0.99 

1.00 

1.14 

0.94 

0.94 

0.91 

0.94 

unst 

18.5 

-6.3 

-45.5 

-8.5 

1.00 

1.00 

1.01 

1.14 

0.93 

0.92 

0.89 

0.95 

unif 

15.5 

-5.0 

-39.1 

-3.1 

1.00 

0.99 

1.00 

1.14 

0.94 

0.94 

0.91 

0.94 

cat.exch 

15.6 

-5.1 

-39.3 

-3.5 

1.00 

0.99 

1.00 

1.15 

0.95 

0.94 

0.91 

0.94 

time.exch 

15.0 

-4.8 

-38.1 

-2.5 

1.00 

0.99 

1.00 

1.14 

0.94 

0.92 

0.91 

0.94 

RC 

15.5 

-5.5 

-41.1 

-5.2 

1.00 

0.99 

1.00 

1.14 

0.94 

0.93 

0.91 

0.93 


indicates correctly specified model and “ ” indicates misspecified model omitting the Z\ predictor 
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structures and the last four ones refer to local odds ratio structures. These seven structures 
are distinguished in terms of the number of parameters being estimated as well as the 
restrictions placed on associations/correlations between the response indicators at different 
time pairs. 

Regarding the GEE for available data, the missing data impact on bias can be noticed 
for all parameters, the higher biases being observed in the hrst intercept, followed by the 
parameter associated with the covariate Z. The impact of the bias on parameter estimates 
was also clearly noticed by the low coverage rates. Still considering the available data it is 
worth noting that the bias for the parameter associated with covariate Z was reduced by 
more than half when the association structure is modeled. All methods being compared 
exhibited bias when their predictive models are incorrectly specihed, although the bias 
for the two intercepts was considerably reduced. For the parameters associated with the 
covariates, the bias of WGEE method was the same magnitude as those provided by 
ordinary GEE, and the performance of DRGEE was slightly superior to MIGEE. 

Independent estimating equations are efficient for the intercept parameters and re¬ 
gression coefficient associated with the baseline covariate X. As expected, the gain in 
efficiency by modeling the association structure occurs only for the time-varying covariate 
Z. Gomparing to the independence structure, the gain in efficiency ranged from 14% on 
average for DRGEE, 23% for multiple imputation and about 17% for WGEE. 

All misspecihed methods presented empirical coverage rates below nominal level al¬ 
though they are somewhat close to the expected value in some cases, particularly for 
MIGEE. 

Table |2] presents simulation results for complete data in addition to correctly specified 
methods. In terms of bias there is no clear distinction between the two approaches, ex¬ 
cept when an unstructured matrix is chosen, which causes an increase in bias, especially 
for WGEE and DRGEE when only the weight model is correctly specified. For the sce¬ 
nario under consideration, the results suggest that the performances of both WGEE and 
DRGEE are slightly superior to multiple imputation, especially for estimates associated 
with covariate X. 

Gompared to independence structure, all other association structures presented smaller 
standard errors for the time-dependent covariate Z. Although they are very similar within 
each method, the largest gain in efficiency, around 23%, is obtained for the MIGEE 
method, followed by WGEE with 19%. For DRGEE estimators it ranged from 14% to 
18%. 

All correctly specified methods showed coverage rates close to the nominal levels for all 
association structures. In terms of empirical bias, relative efficiency and empirical coverage 
there was no clear distinction between the correlation and local odds parametrizations for 
n = 300 subjects. For the scenario under consideration, an exchangeable correlation 
structure or an uniform local odds structure resulted in good marginal mean estimates. 

Simulation results for n = 50,150 and 600 are presented in the Appendix. With 
sample sizes n = 50 or n = 150 subjects, simulation results suggest that the local odds 
parametrization outperforms the correlation coefficient parametrization in terms of bias 
and convergence issues. 
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Table 2: Evaluation criteria for correctly specified models. Results for n = 300 and 
S = 1000 simulations. _ 

Relative Bias Relative Efficiency Empirical Coverage 


Structure 

/3oi 

/3o2 

X 

Z 

/3oi 

/3o2 

X 

z 

/3oi 

P02 

X 


ind 

2.1 

0.1 

-1.1 

-0.5 

1.00 

1.00 

Complete 

1.00 1.00 

0.96 

0.94 

0.95 

0.95 

exch 

0.3 

0.6 

-0.2 

1.3 

0.99 

0.99 

0.99 

1.25 

0.96 

0.95 

0.97 

0.96 

unst 

0.9 

0.0 

-3.4 

-1.7 

1.00 

1.00 

1.01 

1.25 

0.95 

0.94 

0.95 

0.96 

unit 

0.1 

0.7 

0.7 

0.3 

1.00 

1.00 

1.00 

1.27 

0.96 

0.95 

0.94 

0.94 

cat. exch 

0.2 

0.7 

0.6 

-0.2 

1.00 

1.00 

1.00 

1.27 

0.96 

0.95 

0.95 

0.94 

time.exch 

0.2 

0.3 

0.1 

0.3 

1.00 

1.00 

1.00 

1.27 

0.95 

0.95 

0.96 

0.94 

RC 

0.2 

0.3 

-0.8 

-1.9 

1.00 

1.00 

1.00 

1.27 

0.95 

0.94 

0.95 

0.94 

ind 

3.6 

-0.4 

-1.8 

-0.7 

1.00 

WGEE(r+) 

1.00 1.00 1.00 

0.96 

0.95 

0.94 

0.96 

exch 

1.3 

0.3 

0.2 

0.6 

0.99 

0.98 

1.00 

1.19 

0.96 

0.96 

0.96 

0.96 

unst 

6.2 

-1.6 

-6.8 

-4.1 

1.00 

1.00 

1.01 

1.19 

0.95 

0.94 

0.93 

0.96 

unit 

2.4 

-0.1 

-0.3 

-0.4 

0.99 

0.99 

1.00 

1.18 

0.97 

0.96 

0.96 

0.95 

cat. exch 

2.6 

-0.2 

-0.5 

-1.0 

0.99 

0.99 

1.00 

1.18 

0.97 

0.95 

0.96 

0.95 

time.exch 

1.4 

0.3 

1.8 

1.1 

0.99 

0.99 

1.00 

1.19 

0.97 

0.97 

0.96 

0.96 

RC 

2.5 

-0.5 

-2.5 

-1.6 

0.99 

0.99 

1.00 

1.19 

0.96 

0.96 

0.95 

0.94 

ind 

5.3 

-1.3 

-7.1 

-0.4 

1.00 

MIGEE(a:+) 

1.00 1.00 1.00 

0.96 

0.94 

0.94 

0.94 

exch 

3.5 

-0.8 

-5.6 

0.7 

0.99 

0.99 

1.00 

1.22 

0.95 

0.95 

0.96 

0.95 

unst 

4.6 

-1.5 

-9.4 

-3.7 

0.99 

0.99 

1.00 

1.21 

0.95 

0.95 

0.95 

0.94 

unit 

3.9 

-0.8 

-5.4 

0.2 

1.00 

0.99 

1.00 

1.22 

0.96 

0.95 

0.95 

0.94 

cat. exch 

4.2 

-0.9 

-5.7 

-0.7 

1.00 

0.99 

1.00 

1.23 

0.96 

0.95 

0.96 

0.94 

time.exch 

4.2 

-1.0 

-5.2 

-0.2 

1.00 

0.99 

1.00 

1.22 

0.95 

0.96 

0.96 

0.94 

RC 

3.7 

-1.3 

-8.0 

-3.1 

1.00 

1.00 

1.00 

1.23 

0.94 

0.95 

0.95 

0.92 

ind 

1.0 

0.4 

-1.7 

-1.4 

1.00 

DRGEE(x+, 
1.00 1.00 

r+) 

1.00 

0.95 

0.95 

0.93 

0.95 

exch 

-0.8 

0.9 

0.3 

1.4 

0.99 

0.99 

1.00 

1.16 

0.95 

0.95 

0.96 

0.95 

unst 

2.4 

-0.5 

-5.7 

-3.6 

1.00 

0.99 

1.01 

1.16 

0.94 

0.94 

0.94 

0.94 

unit 

0.0 

0.7 

-0.6 

-0.1 

1.00 

0.99 

1.00 

1.16 

0.95 

0.95 

0.95 

0.94 

cat. exch 

0.1 

0.6 

-0.8 

-0.6 

1.00 

0.99 

1.00 

1.17 

0.95 

0.95 

0.95 

0.94 

time.exch 

-0.9 

1.2 

2.0 

1.4 

0.99 

0.99 

1.00 

1.16 

0.95 

0.95 

0.95 

0.94 

RC 

-0.2 

0.3 

-2.3 

-1.8 

1.00 

1.00 

1.00 

1.16 

0.95 

0.94 

0.94 

0.93 

ind 

2.2 

0.0 

-4.7 

-1.5 

1.00 

DRGEE(a;-, 
1.00 1.00 

r+) 

1.00 

0.95 

0.95 

0.95 

0.94 

exch 

0.2 

0.6 

-2.5 

1.3 

0.98 

0.99 

1.00 

1.17 

0.96 

0.96 

0.96 

0.95 

unst 

4.7 

-1.3 

-11.1 

-3.7 

1.00 

1.01 

1.03 

1.18 

0.94 

0.94 

0.94 

0.94 

unit 

1.2 

0.3 

-3.3 

-0.1 

1.00 

0.99 

1.00 

1.17 

0.96 

0.96 

0.97 

0.94 

cat. exch 

1.3 

0.3 

-3.5 

-0.5 

1.00 

0.99 

1.00 

1.18 

0.96 

0.96 

0.97 

0.94 

time.exch 

0.3 

0.7 

-1.2 

1.3 

1.00 

1.00 

1.00 

1.18 

0.96 

0.96 

0.97 

0.94 

RC 

1.2 

-0.1 

-5.5 

-1.6 

1.00 

1.00 

1.00 

1.17 

0.95 

0.96 

0.96 

0.93 

ind 

1.0 

0.2 

-2.5 

-1.8 

1.00 

DRGEE(a;+, 
1.00 1.00 

r-) 

1.00 

0.95 

0.94 

0.93 

0.94 

exch 

-0.4 

0.6 

-1.4 

0.5 

0.99 

0.99 

1.00 

1.15 

0.94 

0.95 

0.96 

0.95 

unst 

2.1 

-0.5 

-5.7 

-5.6 

1.00 

0.99 

1.01 

1.14 

0.93 

0.94 

0.94 

0.95 

unif 

-0.2 

0.7 

-0.6 

0.0 

1.00 

0.99 

1.00 

1.15 

0.95 

0.95 

0.95 

0.94 

cat. exch 

-0.1 

0.6 

-0.8 

-0.5 

1.00 

0.99 

1.00 

1.15 

0.94 

0.95 

0.95 

0.94 

time.exch 

-0.4 

0.8 

0.1 

0.5 

1.00 

0.99 

1.00 

1.14 

0.94 

0.94 

0.95 

0.94 

RC 

0.0 

0.1 

-3.3 

-2.2 

1.00 

1.00 

1.00 

1.14 

0.94 

0.94 

0.94 

0.93 


indicates correctly specified model and “ ” indicates misspecified model omitting the Z\ predictor 
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6 Data Analysis: Functional Classification in Rheumatic 
Mitral Stenosis 

A cohort of 164 patients with rheumatic mitral stenosis who were referred for treatment 
at Hospital das Clinicas of the Federal University of Minas Gerais, Brazil, was selected for 
a mitral valve invervention. Patients were included before intervention and then followed 
up in the outpatient clinic every 4 months according to their clinical status. The first 
three measurements were available for analysis. 

Mitral stenosis is a narrowing of the mitral valve in the heart caused by rheumatic 
disease, which restricts the flow of blood through the valve. The main clinical manifesta¬ 
tion of this disease is shortness of breath, classihed in four categories based on how much 
the patients are limited during physical activity. The response of interest is The New 
York Heart Association (NYHA) Functional Classihcation, that provides a simple way of 
classifying the extent of shortness of breath. Patients with no symptoms and no limitation 
in ordinary physical activity were classihed in class I; slight limitation of physical activity 
in class H; marked limitation of physical activity in class HI; and patients with severe 
limitations resulting in inability to carry on any physical activity without discomfort in 
class IV. Only one patient were classihed into class IV in the followup evaluation and 
hence class IV was combined with class HI in the analysis. Thus, the ordinal response 
was dehned as (1: if class I, 2: if class H, 3 if class HI or class IV). 

Percutaneous mitral valvuloplasty (PMV) is an ehective treatment for stretching the 
stenosed mitral valve. This procedure is carried out by inserting a catheter with a balloon 
at its tip to open the narrowed mitral valve. This procedure causes improvement of the 
functional class in the majority of the patients. A number of patient’s characteristics 
were measured at baseline, such as atrial compliance (Cn: dehned as 1, if < 4, and 0 
otherwise), cardiac rhythm, morphological features of the mitral valve expressed as an 
echocardiographic score, mitral valve area, pressure transmitral gradients, and pulmonary 
artery pressure. Some variables measured after the procedure include the success of the 
procedure to open the mitral valve without complications, long-term event-free survival, 
mitral valve area, pressure transmitral gradients, and pulmonary artery pressure at the 
follow-up appointment. 

This study was characterized by an arbitrary pattern of missing data. The response 
were fully observed for 125 patients, 29 had only the hrst measurement, 1 patient had 
only the hrst data collected, and for 9 patients the second occasion was missing. There 
was no missing data in response at baseline. Some collected variables, such as success 
of the procedure and long term events were responsible for the missingness at followup. 

A baseline covariate of particular interest, atrial compliance, were missing for 54(32.9%) 
of the patients. Among the reasons for not observing such predictor it can be included 
morphological characteristics of the mitral valve and valve calcihcation. Therefore a MAR 
mechanism seems to be a reasonable assumption for this data set. 

Now, the models used for analysis are described. For the ordinal response it was used 
the following proportional odds model 

logit Pr{NYHAitj < j\uit) = /3oj + uj^f3, j = 1,2, t = 1, 2, 3, (14) 

where Uu is the covariate vector at time t, and it is formed by time, atrial compliance, 
echocardiographic score, and success of the procedure. 
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When using WGEE or DRGGE it is necessary to correctly model the missingness 
mechanism in order to obtain consistent estimates of (3. For the missing data process, 
was dehned as the indicator of observing the response NYHAit and Rf was dehned 
as the indicator of observing the baseline atrial compliance. The probability of observing 
the response was modeled as 

logttPriRl = l)=ro+wfr, t = 2,3, (15) 

where the vector included success of the procedure, long term events, and the previous 
indicator of missing data. The model for Rf was specihed as 

logit Pr{Rf = 1) = V’o + (16) 

where the vector wf included the baseline response, echocardiographic score, valve calci¬ 
fication and EGG rhythm. 

For the covariate, the following model was built 

logit PriCn = 1) = "ufT, (17) 

where vj included the baseline right ventricular systolic pressure, mean gradient and 
mitral valve area. A imputation model for NYHA was also specihed that included the 
covariates in ( 1 TT )1 in addition to the baseline right ventricular systolic pressure and the 
response history. 

Here, four methods were used to analyze the data. Results are shown in Table [3l The 
hrst method is the usual GEE method using the available data; the second is the weighted 
method (WGEE) using models ffThll and flTHll for the weights; the third is the multiple 
imputation by chained equation (MIGEE); and the fourth, labeled DRGEE, is the doubly 
robust method using flT^ and flT6|l for the missing data process and flTTl) for the covariate 
models. In order to account for the dependence structure of the repeated measures the 
correlation coefficient and the local odds ratio parametrization were both initially applied 
with the same association structures presented in the simulation. Gonvergence issues were 
observed for unstructured correlation matrix. Motivated by simulation results that sug¬ 
gested better performance of the local odds ratio structures, results are shown only for 
the independence, uniform and category exchangeability local odds ratio structures. The 
independence structure estimates no parameters, while the uniform represents the depen¬ 
dence structure by a single parameter. The category exchangeability assumes constant 
odds ratio among the levels of the response, but different odds ratio between time pairs. 
Thus, the association is explained by three parameters in this structure. Similar results 
were obtained for the other dependence structures they are shown in Appendix. 

The signihcance of time effect indicates the effectiveness of the intervention and im¬ 
provement of functional class over time. Similar coefficients for the two followup occasions 
suggest the major change in functional class occurs right after the valvuloplasty interven¬ 
tion. Total score was non signihcant for all the four methods, although a marginally 
signihcance {p = 0.060) was noted for MIGEE when the category exchangeability odds 
ratio structure was adopted, as an indication that higher scores may be related to cardiac 
insufficiency. All methods provide the same conclusion for effects of the missing covariate 
Gn. It can be noticed that estimates for success of the procedure effect goes from non 
signihcant in the standard GEE to signihcant for all methods, except WGEE, regardless 


16 


Table 3: Results for Rheumatic Mitral Stenosis study under independence, uniform and 
category exchangeability association structures 



Independence 


Uniform 


Cat. Exchangeability 

Parameter 

Est. 

SE 

P 

Est. 

SE 

P 

Est. 

SE 

P 

/^Ol 

-0.800 

0.765 

0.295 

-0.838 

Available 

0.735 

0.254 

-0.768 

0.751 

0.307 

Po2 

0.997 

0.745 

0.181 

0.954 

0.726 

0.189 

1,037 

0.746 

0.164 

Success 

0.733 

0.441 

0.096 

0.662 

0.391 

0.090 

0.692 

0.400 

0.083 

Total Score 

-0.147 

0.104 

0.158 

-0.142 

0.101 

0.161 

-0.156 

0.103 

0.130 

Cn 

0.533 

0.328 

0.104 

0.498 

0.315 

0.113 

0.499 

0.318 

0.116 

time=2 

1,003 

0.413 

0.015 

1,084 

0.389 

0.005 

1,126 

0.390 

0.004 

time=3 

1,216 

0.420 

0.004 

1,214 

0.390 

0.002 

1,198 

0.394 

0.002 

/^Ol 

-0.852 

0.786 

0.278 

-0.883 

WGEE 

0.763 

0.247 

-0.792 

0.776 

0.307 

Po2 

0.895 

0.768 

0.244 

0.882 

0.748 

0.238 

0.974 

0.767 

0.204 

Success 

0.750 

0.449 

0.094 

0.732 

0.405 

0.071 

0.749 

0.416 

0.072 

Total Score 

-0.139 

0.107 

0.192 

-0.136 

0.104 

0.190 

-0.153 

0.106 

0.148 

Cn 

0.524 

0.336 

0.118 

0.464 

0.325 

0.153 

0.472 

0.327 

0.148 

time=2 

1,030 

0.413 

0.013 

1,039 

0.401 

0.010 

1,085 

0.406 

0.007 

time=3 

1,269 

0.419 

0.002 

1,352 

0.404 

0.001 

1,331 

0.410 

0.001 

/^Ol 

-0.856 

0.582 

0.141 

-0.805 

MIGEE 

0.580 

0.165 

-0.788 

0.572 

0.169 

Po2 

0.988 

0.582 

0.090 

1,019 

0.574 

0.076 

1,035 

0.572 

0.070 

Success 

0.938 

0.367 

0.010 

0.922 

0.333 

0.006 

0.980 

0.338 

0.004 

Total Score 

-0.141 

0.080 

0.077 

-0.140 

0.079 

0.075 

-0.147 

0.078 

0.060 

Cn 

0.425 

0.305 

0.163 

0.307 

0.261 

0.239 

0.343 

0.267 

0.199 

time=2 

0.995 

0.334 

0.003 

1,106 

0.308 

0.000 

1,127 

0.330 

0.001 

time=3 

0.969 

0.338 

0.004 

1,068 

0.304 

0.000 

1,079 

0.333 

0.001 

Poi 

-0.793 

0.730 

0.277 

-0.830 

DRGEE 

0.711 

0.243 

-0.726 

0.711 

0.308 

1^02 

0.952 

0.722 

0.187 

0.934 

0.707 

0.186 

1.039 

0.712 

0.145 

Success 

0.825 

0.422 

0.050 

0.806 

0.382 

0.035 

0.841 

0.388 

0.030 

Total Score 

-0.150 

0.099 

0.131 

-0.147 

0.097 

0.132 

-0.165 

0.098 

0.092 

Cn 

0.571 

0.341 

0.094 

0.520 

0.329 

0.113 

0.523 

0.333 

0.116 

time=2 

1.013 

0.499 

0.042 

1.005 

0.512 

0.050 

1.020 

0.525 

0.052 

time=3 

1.135 

0.476 

0.017 

1.185 

0.473 

0.012 

1.145 

0.491 

0.020 
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the adopted association structure. It is interesting to note that the estimated effect are 
increased for all methods, especially for the multiple imputation procedure. So, for exam¬ 
ple, for the uniform structure, the estimated odds of response in class I for subjects with 
success of the procedure compared to patients with suboptimal results were = 2.24 
for DRGEE and = 2.51 for MIGEE. Regarding the association structures it can 
be noticed that the uniform and category exchangeability choices presented some gain in 
efficiency, specially for the success of procedure effect. Similar conclusions were reached 
for the other association structures. 


7 Discussion 


In this paper it was considered a doubly robust estimator for analysis of longitudinal 
data when missingness can occur in a baseline covariate or intermittently in the ordinal 
response. The main objective was to compare the performance of the proposed method in 
terms of bias and efficiency for two different approaches to modeling the cova riance matrix 


of the longitudinal outcome, namely, the correlation coeffic i ent pr oposed by iLipsitz et al 


f I 994 I) and local odds ratio proposed by Touloumis et al. ( 2013 1. Although a complete 
data comparison betw een these two distinct approaches had already appeared elsewhere 
(jNooraee et an . 120141 ) such a comparison with MAR missing data was still in need of 
investigation. 

The covariate design plays an important role in the efficiency of GEE estimators. The 
working independence structure is expected to be efficient for time-s tationary covariates . 
This is no longer true for time-varying covariates and/or missing data ( Lipsitz et all 1994 1. 
Simulation results agreed with the literature, that is, the gain in efficiency by adopting 
more complicated association structures was noticed for the time-varying covariate and it 
somewhat varied through the different methods. For complete data, the standard error for 
the independence was, on average, about 27% larger compared to an uniform structure 
that estimates a single association parameter. The gain in efficiency reached 23% for 
MIGE E. 19% for WGEE an d 18% for DRGEE. 

As Liang fc Zegerl ( 19861 ) pointed out when the assumed correlation is the true corre¬ 
lation, the missing completely at random assumption can be unnecessary. However, this 
is not true when missingness occurs in a co yariate that is MAR given th e response. In 
this case even likelihood methods are biased (iGarpenter fc Kenwardl . 1201311 and the miss¬ 
ingness mechanism must be modeled in order to obtain valid estimates for the regression 
parameters. 

The dependence structures compared here differ in terms of number of parameters 
and restrictions imposed on the correlation/association between levels of ordinal response 
at different time pairs. Although the same definition can be used for an independence, 
exchangeable and unstructured associa tion matrix, this does not imply that identical as¬ 
sociations are fit ( Nooraee et all 2014 1. Under the independence working assumption all 
off-diagonal blocks of the covariance matrix are constant and equal to zero. Exchange¬ 
ability over time indicates that the association between Y^j and Yn/ji is independent of 
time, but it depends on the levels j and j'. For unstructured associations there are no 
restrictions implied. With moderate to large number of subjects (at least 300 subjects) 
the simulation results did not suggest a very clear distinction in terms of bias, relative 
efficiency and empirical coverage between the correlation or local odds ratio parametriza- 
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tions, although some bias was observed for unstructured matrices. On the other side, 
for small sample sizes {n = 50 or 150) simulation results did suggest that the local odds 
ratio outperforms the correlation coefficient parametrization in terms of bias and conver¬ 
gence issues. It seems that local odds structures works fine for small samples sizes and 
the correlation coefficient needs relatively more subjects to achieve the same reduction in 
bias. 

It is noticeable two important differences between the local odds ratio and the cor¬ 
relation coefficient parametrizations. Firstly, the local odds ratios do es not depend on 


the m arginal specification, that is, (3 and cx are variation independent fiTonlonmis et al 
2 OI 3 I) . Unlike correlations, the local odds estimates for the association structure does not 


depend on the values of the observed covariates and thus these estimates do not need to 
be obtained at each step of the modified Fisher Scoring algorithm. As a consequence, 
the iterative procedure converges faster. Secondly, the correlation estimates are obtained 
through a moment based estimator, thus some work was necessary to ensure their validity 
in the presence of data that is MAR. The estimates of local odds are based on maximum 
likelihood methods, assuming independent Poisson sampling for the observed counts in the 
marginalized contingency tables. As it considers only the observed responses the resulting 
estimates are valid under a MAR mechanism. 

In some situations (small samples, complex patterns of missing data), it is possible 
that the algorithm will not conv erge because the es timated correlation matrix is not 
guaranteed to be positive definite (iLipsitz et Il994l ). For small sample sizes there may 
not be enough data to estimate both the regression parameters and a correlation matrix 
that is highly unstructured. Simulation results for unstructured correlation matrices and 
small sample sizes presented very low convergence rates (in addition to increased bias) and 
thus an exchangeable correlation matrix is preferred. There were no serious convergence 
issues with the local odds ratio parametrization, except for an unstructured matrix applied 
to n = 50. 

Although the latent vectors have been generated from an exchangeable correlation, the 
correlation coefficient betwe en binary variables Y ni and Uf'y also depends on the linear 


predictor in times t and t' f Nooraee et al . 2014ll . That is, the correlation between the 


latent vectors was exchangeable, but the correlation between the binary variables was 
not exchangeable due to the dependent mean in the marginal model. Nevertheless, it is 
expected the GEE method remain valid even assu ming a misspecified corre l ation matrix 
provided the marginal mean is correctly specified f Molenberghs fc Kenward . 2 OIC 1 I 1 . 

Overall, the doubly robust method performed well for all working association structures 
under comparison except when applied to a very small number of subjects (i.e., n = 50) 
and an unstructured matrix was adopted. The coverage probabilities were relatively close 
to the nominal level of 95%. 

The doubly robust estimator considered here is restricted to missingness in the response 
and a baseline covariate. A natural extension is to consider an intermittently missing time- 
varying covariate and/or allow multiple missing covariates. In the proposed approach the 
marginal means were modeled by cumulative logits. This implies a proportional odds 
model that in some cases may not be valid. Another possible extension of the proposed 
model is, therefore, to a l low n on-proportional odds for a subset of the explanatory variables 
( Peterson fc Harrell Jr . 1990 b 
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A Appendix 

A.l Asymptotic Variance 

To state the asymptotic properties of /3, let 
S'ij(/3, 7 ) be the individual’s contribution to the estimating equations for (3, 

‘S' 2 i(V’) be the individual’s contribution to the estimating equations for -i/j, and 
be the individual’s contribution to the estimating equations for 7 . 

Dehne = E /i2(/3,V’,7) = E 

/i3(/3, V^,7) = E{dSM^P,7)/d-f^}, = E{dS2^W/d^P^}, hh) = E{dSs,i^)/d-f^}, 

and Qi(/3,V’,7) = - Ii2{/3,'ip,j)l2\'ip)S2i{i’) - Ii3{/3,'ip,-f)l3\j)S3ih). 


Theorem 1 If either the missing data model or the covariate model is correctly specified, 
then 


nV2(/3 - /3„) iV(0, r-i(/3o, 7o)S {r-'(/3o, ^Po^ 7o)} ), 


(18) 


where is the true value of (3, ipQ and Jq are the probability limits of xp and 7 , and 
Yi = E |Qj(/3Q, xpQ, 7o)Qj (/3o, V’oj 7o)} ■ 

Inferences for [3 follows by replacing the unknown quantitie s in flT^ bv i ts consistent 
estimators. We make use of “generalized information equality” fjPiercel . 1198211 that 
E{dSu{f3,xP,l)/dxP^] = -E{Su{f3,xP,-t)Sl{xP)], and ^^ 

E ■[9S'ij(/3, xp, 7)/97^} = —E |5'ii(/3, xp, 7 ) 5 ^ 34 ( 7 )}. Similarly ( Robins et all 1995 1. 

E {dS 2 iixP)/dxP^} = -Var{S 2 iixP)}, and E {dS 3 iij)/dj^} = -Var {Ss^ij)}. 


The matrix T is replaced by F = n ^ \ ^Eu^O) jdfB and S by S = n ^ Y^^=i \ QiQ 




-1 , 


Q, = Su{e)-I,2{e)I^ ixP)S2^ixP)-h3i9)l3 (7)-53*(7), luiO) = Er=i dSuie)/dxP 


1,3(6) = n-i Er=i {dSu(e)/dY]. W) = n-^ ELi ^ dS2.(xP)/dxP 
h(l) = Er=i {dS3i{ ^)ld-i'^]. 


The proof is similar to I Chen fc Zhoul (1201111 and is omitted here. 
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A.2 Additional Simulation Results 


Table 4: Evaluation criteria for misspecified models. Results for n = 50 and S = 1000 
simulations._ 




Relative Bias 


Relative Efficiency 


Empirical Coverage 


Structure 

/3oi 

/3o2 

X 

Z 

/3oi 

/3o2 

X 

Z 

/3oi 

/3o2 

X 

z 

ind 

101.9 

-25.2 

-10.0 

48.4 

1.00 

Available 
1.00 1.00 

1.00 

0.85 

0.85 

1.00 

0.95 

exch 

96.2 

-22.6 

-36.0 

19.5 

1.02 

1.03 

1.03 

1.15 

0.87 

0.88 

1.00 

0.95 

unst 

91.0 

-21.8 

-44.1 

-1.3 

1.08 

1.10 

1.10 

1.17 

0.86 

0.90 

1.00 

0.96 

unif 

80.2 

-18.6 

-14.8 

26.5 

1.01 

0.99 

1.00 

1.19 

0.91 

0.90 

1.00 

0.94 

cat. exch 

80.9 

-19.0 

-16.8 

19.5 

1.01 

1.00 

1.01 

1.23 

0.90 

0.89 

1.00 

0.94 

time.exch 

82.8 

-18.2 

-16.9 

25.1 

1.00 

0.99 

1.01 

1.20 

0.88 

0.90 

1.00 

0.93 

RC 

83.6 

-20.0 

-19.2 

9.1 

1.01 

1.00 

1.01 

1.26 

0.95 

0.93 

1.00 

0.93 

ind 

24.5 

-0.7 

-13.5 

30.3 

1.00 

WGEE(r-) 
1.00 1.00 

1.00 

0.94 

0.93 

1.00 

0.96 

exch 

36.3 

-5.2 

-29.7 

14.2 

1.01 

1.02 

1.03 

1.14 

0.94 

0.94 

1.00 

0.95 

unst 

43.3 

-7.0 

-44.7 

-3.0 

1.07 

1.11 

1.09 

1.17 

0.93 

0.93 

1.00 

0.96 

unif 

17.9 

0.6 

-12.6 

24.8 

1.00 

0.99 

0.99 

1.14 

0.96 

0.94 

1.00 

0.96 

cat.exch 

19.2 

0.0 

-11.6 

19.0 

1.00 

0.99 

0.99 

1.15 

0.95 

0.93 

1.00 

0.95 

time.exch 

18.1 

0.6 

-12.6 

23.6 

0.98 

0.98 

0.99 

1.15 

0.96 

0.95 

1.00 

0.94 

RC 

22.3 

-1.4 

-12.7 

9.9 

0.99 

0.99 

0.99 

1.17 

0.97 

0.95 

1.00 

0.95 

ind 

25.6 

-5.0 

-34.7 

-2.2 

1.00 

MIGEE(a:-) 
1.00 1.00 

1.00 

0.96 

0.95 

1.00 

0.96 

exch 

35.4 

-7.8 

-53.5 

-8.3 

1.00 

1.01 

1.00 

1.18 

0.95 

0.97 

1.00 

0.94 

unst 

35.2 

-6.7 

-59.9 

-28.8 

1.03 

1.06 

1.05 

1.16 

0.96 

0.97 

1.00 

0.95 

unif 

23.0 

-5.1 

-36.4 

2.4 

1.00 

0.99 

1.00 

1.19 

0.97 

0.96 

1.00 

0.94 

cat. exch 

24.8 

-5.8 

-35.3 

-3.5 

0.99 

0.99 

1.00 

1.24 

0.96 

0.96 

1.00 

0.95 

time.exch 

22.9 

-4.0 

-36.0 

-1.2 

0.98 

0.99 

0.99 

1.21 

0.97 

0.96 

1.00 

0.92 

RC 

25.5 

-5.5 

-35.8 

-16.3 

1.00 

0.99 

1.00 

1.28 

0.95 

0.94 

1.00 

0.93 

ind 

21.1 

-1.4 

-19.0 

-4.0 

DRGEE(; 
1.00 1.00 

x ~, r~ 

1.00 

■) 

1.00 

0.94 

0.93 

1.00 

0.95 

exch 

33.8 

-6.0 

-35.9 

-14.4 

1.02 

1.03 

1.04 

1.18 

0.93 

0.94 

1.00 

0.92 

unst 

38.2 

-7.8 

-55.3 

-34.2 

1.08 

1.13 

1.13 

1.20 

0.93 

0.93 

1.00 

0.93 

unif 

15.4 

-0.1 

-18.2 

-1.0 

1.00 

0.98 

0.99 

1.15 

0.96 

0.94 

1.00 

0.94 

cat. exch 

16.3 

-0.7 

-18.9 

-4.8 

1.00 

0.98 

0.99 

1.17 

0.95 

0.93 

1.00 

0.94 

time.exch 

15.5 

-0.1 

-17.1 

-2.5 

0.99 

0.98 

1.00 

1.17 

0.96 

0.94 

1.00 

0.92 

RC 

20.1 

-1.7 

-20.1 

-15.2 

0.99 

0.98 

0.99 

1.21 

0.95 

0.94 

1.00 

0.91 


indicates correctly specified model and “ ” indicates misspecified model omitting the Zi predictor 
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Table 5: Evaluation criteria for correctly specified models. Results for n = 50 and S 
1000 simulations._ 


Relative Bias Relative Efficiency Empirical Coverage 


Structure 

/3oi 

P02 

X 

Z 

/3oi 

/3o2 

X 

Z 

/3oi 

/3o2 

X 


ind 

3.4 

4.4 

3.3 

3.5 

1.00 

Complete 
1.00 1.00 

1.00 

0.95 

0.94 

0.99 

0.96 

exch 

13.6 

1.3 

-12.4 

-8.0 

0.99 

1.01 

1.01 

1.23 

0.94 

0.96 

1.00 

0.94 

unst 

17.0 

2.3 

-25.5 

-18.9 

1.03 

1.06 

1.05 

1.21 

0.95 

0.96 

0.98 

0.94 

unit 

0.1 

4.6 

4.0 

2.8 

1.00 

1.00 

1.00 

1.25 

0.95 

0.93 

1.00 

0.95 

cat. exch 

0.9 

4.2 

3.6 

-1.8 

1.00 

1.00 

1.00 

1.28 

0.95 

0.93 

0.99 

0.95 

time.exch 

3.3 

4.8 

-2.6 

1.6 

0.99 

0.99 

1.00 

1.26 

0.95 

0.95 

0.99 

0.93 

RC 

5.1 

3.0 

-0.6 

-8.1 

1.00 

1.00 

1.00 

1.31 

0.95 

0.93 

1.00 

0.93 

ind 

14.8 

3.2 

1.8 

9.1 

1.00 

WGEE(r+) 
1.00 1.00 

1.00 

0.95 

0.94 

1.00 

0.97 

exch 

27.7 

-1.7 

-16.3 

-1.6 

1.01 

1.03 

1.03 

1.16 

0.94 

0.94 

1.00 

0.95 

unst 

37.5 

-4.7 

-29.3 

-13.5 

1.08 

1.12 

1.10 

1.18 

0.93 

0.93 

1.00 

0.97 

unit 

8.1 

5.1 

4.7 

10.1 

1.00 

0.98 

0.99 

1.15 

0.96 

0.93 

1.00 

0.95 

cat. exch 

10.0 

4.4 

4.4 

4.9 

1.00 

0.98 

1.00 

1.16 

0.96 

0.94 

1.00 

0.95 

time.exch 

6.9 

5.7 

8.6 

5.7 

0.97 

0.97 

0.99 

1.16 

0.96 

0.95 

1.00 

0.95 

RC 

12.4 

2.9 

0.8 

-5.9 

0.98 

0.98 

0.99 

1.19 

0.96 

0.94 

1.00 

0.94 

ind 

12.1 

0.6 

-5.1 

6.5 

1.00 

MIGEE(x+) 
1.00 1.00 

1.00 

0.96 

0.95 

1.00 

0.96 

exch 

23.7 

-2.9 

-20.3 

-7.3 

0.99 

1.01 

1.01 

1.18 

0.95 

0.96 

1.00 

0.94 

unst 

26.8 

-3.3 

-32.4 

-27.0 

1.02 

1.06 

1.04 

1.15 

0.95 

0.97 

0.99 

0.95 

unit 

7.5 

0.9 

0.0 

5.1 

1.00 

1.00 

1.00 

1.19 

0.96 

0.94 

1.00 

0.94 

cat. exch 

9.1 

0.3 

0.5 

-0.6 

1.00 

1.00 

1.00 

1.23 

0.95 

0.94 

1.00 

0.94 

time.exch 

7.8 

1.6 

-2.1 

2.1 

0.99 

0.99 

1.00 

1.21 

0.96 

0.95 

1.00 

0.93 

RC 

13.4 

-0.1 

-4.9 

-14.0 

1.00 

0.99 

1.00 

1.28 

0.95 

0.94 

1.00 

0.92 

ind 

7.0 

4.1 

0.1 

2.4 

1.00 

DRGEE(x+,r+) 

1.00 1.00 1.00 

0.94 

0.94 

1.00 

0.96 

exch 

19.3 

-0.5 

-13.9 

-10.0 

1.01 

1.03 

1.04 

1.19 

0.94 

0.95 

1.00 

0.93 

unst 

25.0 

-2.5 

-28.3 

-23.3 

1.08 

1.13 

1.13 

1.21 

0.92 

0.94 

1.00 

0.94 

unit 

3.9 

5.5 

3.7 

4.5 

0.99 

0.98 

0.99 

1.16 

0.95 

0.93 

1.00 

0.94 

cat. exch 

4.4 

4.8 

3.9 

-0.1 

0.99 

0.98 

0.99 

1.18 

0.94 

0.93 

1.00 

0.94 

time.exch 

2.7 

5.7 

10.5 

0.8 

0.97 

0.97 

0.99 

1.16 

0.96 

0.95 

1.00 

0.93 

RC 

9.0 

3.2 

1.0 

-10.7 

0.97 

0.97 

0.98 

1.20 

0.95 

0.94 

1.00 

0.91 

ind 

11.6 

3.3 

-2.4 

2.3 

1.00 

DRGEE(a;-, r 
1.00 1.00 

1.00 

0.94 

0.94 

1.00 

0.96 

exch 

25.8 

-2.0 

-22.2 

-9.0 

1.02 

1.06 

1.06 

1.22 

0.93 

0.95 

1.00 

0.93 

unst 

31.8 

-5.0 

-39.5 

-26.1 

1.09 

1.18 

1.18 

1.26 

0.92 

0.93 

1.00 

0.93 

unit 

5.4 

5.1 

-1.3 

3.2 

1.00 

0.98 

1.01 

1.17 

0.95 

0.94 

1.00 

0.95 

cat. exch 

7.3 

4.5 

-2.7 

-1.3 

1.00 

0.98 

1.01 

1.20 

0.95 

0.94 

1.00 

0.94 

time.exch 

5.2 

5.2 

1.6 

1.1 

0.98 

0.97 

0.99 

1.18 

0.96 

0.95 

1.00 

0.94 

RC 

11.2 

3.1 

-4.7 

-12.5 

0.98 

0.97 

1.00 

1.22 

0.94 

0.94 

1.00 

0.92 

ind 

8.9 

3.9 

-2.0 

2.0 

1.00 

DRGEE(a;+,r- 
1.00 1.00 

-) 

1.00 

0.93 

0.93 

1.00 

0.96 

exch 

18.3 

-0.5 

-14.8 

-12.9 

1.01 

1.03 

1.03 

1.18 

0.94 

0.95 

1.00 

0.92 

unst 

23.7 

-2.3 

-28.5 

-30.3 

1.06 

1.11 

1.11 

1.20 

0.92 

0.94 

1.00 

0.93 

unit 

4.7 

4.9 

3.7 

2.6 

0.99 

0.98 

0.99 

1.15 

0.95 

0.93 

1.00 

0.95 

cat. exch 

5.7 

4.2 

2.4 

-2.4 

0.99 

0.98 

0.99 

1.17 

0.95 

0.92 

1.00 

0.94 

time.exch 

3.3 

5.3 

4.4 

1.3 

0.98 

0.98 

1.00 

1.16 

0.96 

0.94 

1.00 

0.92 

RC 

10.5 

3.4 

-2.2 

-11.9 

0.99 

0.98 

0.99 

1.19 

0.88 

0.89 

1.00 

0.92 


indicates correctly specified model and “ ” indicates misspecified model omitting the Zi predictor 
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Table 6: Evaluation criteria for misspecified models. Results for n = 150 and S = 1000 
simulations._ 




Relative Bias 


Relative Efficiency 


Empirical Coverage 


Structure 

/3oi 

/3o2 


Z 

/3oi 

/3o2 


Z 

Poi 

/3o2 

X 

Z 

ind 

117.0 

-29.9 

-20.9 

60.5 

1.00 

Available 

1.00 1.00 

1.00 

0.46 

0.66 

0.94 

0.84 

exch 

96.4 

-25.5 

-32.2 

30.2 

1.02 

1.01 

1.02 

1.20 

0.58 

0.7 

0.92 

0.92 

unst 

97.4 

-26.3 

-40.0 

17.4 

1.04 

1.03 

1.04 

1.21 

0.59 

0.71 

0.94 

0.94 

unif 

92.5 

-22.6 

-25.9 

26.8 

1.02 

1.00 

1.01 

1.23 

0.63 

0.78 

0.93 

0.91 

cat. exch 

92.6 

-22.7 

-26.4 

24.8 

1.02 

1.00 

1.01 

1.24 

0.64 

0.77 

0.93 

0.92 

time.exch 

93.6 

-22.8 

-25.6 

30.0 

1.02 

1.00 

1.01 

1.23 

0.61 

0.76 

0.94 

0.91 

RC 

93.2 

-22.4 

-25.1 

24.1 

1.02 

1.00 

1.01 

1.25 

0.62 

0.78 

0.94 

0.92 

ind 

19.0 

-3.4 

-22.9 

37.5 

1.00 

WGEE(r-) 
1.00 1.00 

1.00 

0.94 

0.94 

0.94 

0.9 

exch 

21.5 

-5.5 

-32.9 

22.2 

1.00 

1.00 

1.00 

1.16 

0.94 

0.93 

0.94 

0.94 

unst 

27.3 

-7.4 

-44.2 

10.0 

1.02 

1.02 

1.02 

1.17 

0.94 

0.94 

0.95 

0.94 

unif 

18.7 

-3.4 

-24.0 

24.7 

1.00 

1.00 

1.00 

1.16 

0.94 

0.94 

0.94 

0.93 

cat. exch 

19.1 

-3.6 

-24.1 

23.0 

1.00 

1.00 

1.00 

1.17 

0.94 

0.94 

0.94 

0.93 

time.exch 

19.0 

-4.0 

-25.3 

26.4 

1.00 

1.00 

1.00 

1.16 

0.94 

0.94 

0.94 

0.92 

RC 

20.3 

-4.0 

-24.8 

20.9 

1.00 

1.00 

1.00 

1.16 

0.94 

0.96 

0.94 

0.94 

ind 

21.2 

-5.9 

-46.4 

-7.4 

1.00 

MIGEE(x-) 
1.00 1.00 

1.00 

0.94 

0.95 

0.94 

0.95 

exch 

22.7 

-7.5 

-53.4 

-2.2 

0.98 

0.99 

0.99 

1.21 

0.95 

0.95 

0.95 

0.96 

unst 

24.9 

-8.1 

-57.2 

-10.0 

0.99 

1.00 

1.00 

1.21 

0.95 

0.96 

0.96 

0.94 

unif 

20.8 

-5.9 

-45.9 

-1.8 

1.00 

1.00 

1.00 

1.22 

0.94 

0.95 

0.94 

0.93 

cat. exch 

21.3 

-6.1 

-46.4 

-3.7 

1.00 

1.00 

1.00 

1.24 

0.94 

0.95 

0.94 

0.93 

time.exch 

20.1 

-6.2 

-44.8 

0.4 

1.00 

1.00 

1.00 

1.22 

0.95 

0.96 

0.96 

0.94 

RC 

20.8 

-6.1 

-45.3 

-3.8 

1.00 

1.00 

1.00 

1.25 

0.95 

0.96 

0.96 

0.94 

ind 

16.6 

-3.5 

-32.9 

-5.4 

DRGEE(; 
1.00 1.00 

x ~, r~ 
1.00 

■) 

1.00 

0.93 

0.94 

0.94 

0.94 

exch 

19.4 

-5.5 

-42.8 

-4.9 

1.00 

1.00 

1.01 

1.13 

0.93 

0.93 

0.93 

0.95 

unst 

24.7 

-7.4 

-54.6 

-16.1 

1.01 

1.02 

1.03 

1.16 

0.93 

0.93 

0.93 

0.93 

unif 

16.0 

-3.3 

-31.9 

-1.7 

1.00 

1.00 

0.99 

1.14 

0.93 

0.94 

0.95 

0.92 

cat.exch 

16.4 

-3.5 

-32.2 

-2.9 

1.00 

1.00 

0.99 

1.16 

0.94 

0.94 

0.94 

0.92 

time.exch 

16.0 

-3.8 

-34.2 

0.6 

1.00 

1.00 

1.00 

1.16 

0.93 

0.94 

0.93 

0.95 

RC 

17.0 

-4.0 

-35.3 

-5.7 

1.00 

1.00 

1.00 

1.17 

0.94 

0.95 

0.94 

0.94 


indicates correctly specified model and “ ” indicates misspecified model omitting the Z\ predictor 
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Table 7: Evaluation criteria for correctly specified models. Results for n = 150 and 
S = 1000 simulations. _ 

Relative Bias Relative Efficiency Empirical Coverage 


Structure 

Poi 

P02 

X 

Z 

/3oi 

/3o2 

X 

Z 

/3oi 

/3o2 

X 


ind 

0.2 

1.9 

2.9 

1.0 

1.00 

Complete 
1.00 1.00 

1.00 

0.94 

0.94 

0.95 

0.95 

exch 

0.9 

0.8 

-3.8 

0.7 

0.99 

1.00 

1.00 

1.26 

0.93 

0.93 

0.95 

0.97 

unst 

3.4 

0.4 

-8.5 

-6.1 

1.01 

1.01 

1.01 

1.26 

0.95 

0.95 

0.97 

0.95 

unit 

-0.1 

2.0 

3.2 

1.6 

1.00 

1.00 

1.00 

1.26 

0.94 

0.94 

0.94 

0.95 

cat. exch 

0.0 

1.9 

2.8 

0.3 

1.00 

1.00 

1.00 

1.27 

0.94 

0.94 

0.94 

0.94 

time.exch 

0.3 

1.6 

3.3 

2.9 

1.00 

1.00 

1.00 

1.27 

0.95 

0.96 

0.95 

0.94 

RC 

0.7 

1.8 

1.3 

-0.1 

1.01 

1.01 

1.00 

1.28 

0.95 

0.95 

0.95 

0.96 

ind 

3.3 

1.5 

3.9 

6.0 

1.00 

WGEE(r+) 
1.00 1.00 

1.00 

0.96 

0.95 

0.96 

0.95 

exch 

5.7 

-0.5 

-7.1 

1.2 

0.99 

0.99 

1.00 

1.16 

0.94 

0.94 

0.95 

0.96 

unst 

13.3 

-2.8 

-18.2 

-7.9 

1.02 

1.02 

1.02 

1.18 

0.96 

0.97 

0.97 

0.95 

unit 

2.6 

1.6 

4.0 

3.7 

0.99 

0.99 

0.99 

1.17 

0.95 

0.95 

0.96 

0.94 

cat. exch 

3.0 

1.4 

3.9 

2.5 

0.99 

0.99 

0.99 

1.18 

0.95 

0.96 

0.96 

0.94 

time.exch 

2.1 

1.3 

2.1 

4.4 

0.98 

0.98 

0.99 

1.17 

0.95 

0.95 

0.95 

0.94 

RC 

4.4 

0.9 

1.5 

-0.3 

1.00 

0.99 

0.99 

1.17 

0.95 

0.96 

0.96 

0.95 

ind 

5.4 

0.1 

-3.3 

0.9 

1.00 

MIGEE(x+) 
1.00 1.00 

1.00 

0.94 

0.95 

0.95 

0.95 

exch 

7.7 

-1.8 

-14.5 

0.3 

0.99 

0.99 

0.99 

1.20 

0.94 

0.94 

0.96 

0.96 

unst 

9.6 

-2.6 

-18.9 

-7.5 

1.00 

1.00 

1.00 

1.20 

0.96 

0.96 

0.96 

0.95 

unit 

5.0 

0.1 

-3.0 

1.4 

1.00 

1.00 

1.00 

1.21 

0.94 

0.95 

0.95 

0.93 

cat. exch 

5.6 

-0.1 

-3.8 

-0.6 

1.00 

1.00 

1.00 

1.23 

0.94 

0.95 

0.95 

0.93 

time.exch 

4.7 

-0.3 

-4.1 

3.3 

1.00 

1.00 

1.00 

1.22 

0.95 

0.95 

0.95 

0.94 

RC 

5.3 

0.0 

-3.7 

-0.8 

1.00 

1.00 

1.00 

1.24 

0.95 

0.95 

0.96 

0.94 

ind 

0.3 

2.5 

5.0 

1.9 

1.00 

DRGEE(x+,r+) 

1.00 1.00 1.00 

0.95 

0.94 

0.95 

0.94 

exch 

3.0 

0.3 

-5.7 

-0.1 

0.98 

0.99 

1.00 

1.14 

0.93 

0.94 

0.95 

0.94 

unst 

8.9 

-1.7 

-16.5 

-9.6 

1.00 

1.01 

1.03 

1.16 

0.94 

0.95 

0.97 

0.93 

unit 

-0.2 

2.5 

5.5 

1.9 

1.00 

1.00 

1.00 

1.15 

0.95 

0.94 

0.95 

0.93 

cat. exch 

0.2 

2.3 

5.1 

0.7 

1.00 

1.00 

1.00 

1.16 

0.94 

0.94 

0.95 

0.92 

time.exch 

-0.6 

2.2 

3.6 

4.3 

0.98 

0.98 

0.99 

1.15 

0.94 

0.95 

0.95 

0.94 

RC 

1.7 

1.8 

2.0 

-1.5 

0.99 

1.00 

0.99 

1.17 

0.95 

0.95 

0.96 

0.94 

ind 

1.8 

1.9 

0.7 

1.9 

1.00 

DRGEE(a;-, r 
1.00 1.00 

+■) 

1.00 

0.95 

0.95 

0.97 

0.94 

exch 

4.6 

-0.2 

-10.2 

-0.7 

0.99 

1.00 

1.01 

1.17 

0.94 

0.94 

0.96 

0.95 

unst 

11.4 

-2.5 

-22.0 

-10.3 

1.02 

1.04 

1.05 

1.19 

0.95 

0.96 

0.96 

0.93 

unit 

1.3 

2.1 

0.8 

1.8 

1.00 

0.99 

1.00 

1.18 

0.95 

0.96 

0.97 

0.94 

cat. exch 

1.7 

1.9 

0.2 

0.6 

1.00 

0.99 

1.01 

1.19 

0.95 

0.96 

0.97 

0.93 

time.exch 

1.0 

1.8 

0.3 

4.4 

0.99 

0.99 

0.99 

1.18 

0.94 

0.95 

0.97 

0.94 

RC 

3.1 

1.4 

-1.6 

-1.8 

1.00 

1.00 

1.00 

1.20 

0.95 

0.95 

0.97 

0.94 

ind 

1.2 

2.1 

3.0 

1.7 

1.00 

DRGEE(a;+,r- 
1.00 1.00 

-) 

1.00 

0.95 

0.94 

0.95 

0.94 

exch 

3.0 

0.3 

-6.6 

-1.6 

1.00 

1.00 

1.01 

1.13 

0.93 

0.94 

0.94 

0.95 

unst 

8.1 

-1.6 

-16.0 

-13.3 

1.00 

1.01 

1.03 

1.16 

0.94 

0.95 

0.97 

0.93 

unit 

0.7 

2.2 

3.4 

1.6 

1.00 

1.00 

1.00 

1.14 

0.94 

0.94 

0.95 

0.93 

cat. exch 

1.1 

2.0 

2.6 

0.4 

1.00 

1.00 

1.00 

1.15 

0.94 

0.94 

0.94 

0.93 

time.exch 

0.3 

1.8 

0.6 

3.7 

1.00 

1.00 

1.00 

1.15 

0.93 

0.94 

0.95 

0.95 

RC 

1.5 

1.6 

2.0 

-2.5 

1.00 

1.00 

1.00 

1.17 

0.94 

0.95 

0.96 

0.94 


indicates correctly specified model and “ ” indicates misspecified model omitting the Z\ predictor 
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Table 8: Evaluation criteria for misspecified models. Results for n = 600 and S = 1000 
simulations._ 




Relative Bias 


Relative Efficiency 


Empirical Coverage 


Structure 

/3oi 

P02 

X 

Z 

Poi 

/3o2 

X 

Z 

Poi 

P02 

X 

Z 

ind 

119.9 

-31.8 

-22.5 

61.8 

1.00 

Available 

1.00 1.00 

1.00 

0.01 

0.08 

0.92 

0.50 

exch 

95.2 

-24.9 

-25.1 

30.8 

1.02 

1.01 

1.01 

1.21 

0.08 

0.27 

0.89 

0.76 

unst 

94.3 

-26.2 

-27.9 

29.8 

1.03 

1.02 

1.02 

1.21 

0.09 

0.20 

0.88 

0.79 

unif 

93.5 

-24.5 

-28.5 

25.3 

1.02 

1.00 

1.01 

1.24 

0.08 

0.26 

0.90 

0.82 

cat. exch 

93.5 

-24.5 

-28.5 

24.7 

1.02 

1.00 

1.01 

1.24 

0.08 

0.26 

0.90 

0.83 

time.exch 

94.9 

-24.7 

-29.4 

25.3 

1.02 

1.00 

1.01 

1.25 

0.07 

0.28 

0.89 

0.83 

RC 

93.4 

-24.0 

-26.9 

26.0 

1.02 

1.00 

1.01 

1.24 

0.08 

0.31 

0.89 

0.82 

ind 

19.1 

-5.6 

-28.1 

34.8 

1.00 

WGEE(r-) 
1.00 1.00 

1.00 

0.92 

0.92 

0.91 

0.82 

exch 

18.5 

-5.6 

-28.5 

21.7 

1.00 

0.98 

0.99 

1.18 

0.92 

0.92 

0.90 

0.89 

unst 

19.0 

-6.0 

-31.2 

20.3 

1.00 

0.99 

0.99 

1.18 

0.92 

0.93 

0.90 

0.88 

unif 

19.2 

-6.0 

-30.3 

20.8 

1.00 

1.00 

0.99 

1.18 

0.92 

0.92 

0.91 

0.90 

cat. exch 

19.3 

-6.0 

-30.3 

20.4 

1.00 

1.00 

0.99 

1.18 

0.92 

0.92 

0.91 

0.90 

time.exch 

19.8 

-6.0 

-30.9 

20.8 

1.00 

1.00 

0.99 

1.18 

0.92 

0.91 

0.91 

0.90 

RC 

18.7 

-5.2 

-27.9 

22.1 

1.00 

0.99 

0.99 

1.18 

0.92 

0.92 

0.89 

0.88 

ind 

19.9 

-7.2 

-47.9 

-8.1 

1.00 

MIGEE(x-) 
1.00 1.00 

1.00 

0.92 

0.91 

0.86 

0.92 

exch 

20.5 

-7.0 

-46.7 

-2.6 

1.00 

1.00 

1.00 

1.21 

0.91 

0.90 

0.87 

0.94 

unst 

20.5 

-7.4 

-48.8 

-3.2 

1.00 

1.00 

0.99 

1.22 

0.90 

0.90 

0.84 

0.96 

unif 

19.8 

-7.3 

-47.8 

-3.7 

1.00 

1.00 

0.99 

1.22 

0.91 

0.90 

0.86 

0.95 

cat. exch 

19.8 

-7.3 

-47.8 

-4.2 

1.00 

1.00 

0.99 

1.22 

0.92 

0.90 

0.86 

0.94 

time.exch 

20.3 

-7.4 

-48.5 

-3.7 

1.00 

1.01 

1.00 

1.22 

0.91 

0.89 

0.83 

0.94 

RC 

19.7 

-6.9 

-46.6 

-3.3 

1.00 

1.00 

0.99 

1.22 

0.90 

0.90 

0.85 

0.95 

ind 

16.1 

-5.8 

-41.1 

-6.7 

1.00 

DRGEE( 

1.00 

x ~, r~ 
1.00 

■) 

1.00 

0.92 

0.91 

0.86 

0.92 

exch 

15.8 

-5.5 

-39.7 

-2.9 

0.99 

0.98 

0.99 

1.14 

0.93 

0.90 

0.86 

0.94 

unst 

16.4 

-5.9 

-42.1 

-3.9 

0.99 

0.98 

0.99 

1.14 

0.92 

0.92 

0.86 

0.96 

unif 

15.7 

-5.8 

-41.1 

-3.4 

0.99 

0.99 

0.99 

1.15 

0.92 

0.91 

0.86 

0.95 

cat.exch 

15.8 

-5.9 

-41.1 

-3.8 

0.99 

0.99 

0.99 

1.15 

0.92 

0.91 

0.86 

0.95 

time.exch 

16.3 

-5.9 

-41.5 

-3.5 

0.99 

0.99 

0.99 

1.15 

0.92 

0.90 

0.87 

0.94 

RC 

15.2 

-5.2 

-38.6 

-2.0 

1.00 

0.99 

0.99 

1.14 

0.92 

0.91 

0.86 

0.95 


indicates correctly specified model and “ ” indicates misspecified model omitting the Z\ predictor 
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Table 9: Evaluation criteria for correctly specified models. Results for n = 600 and 
S = 1000 simulations. _ 

Relative Bias Relative Efficiency Empirical Coverage 


Structure 

/3oi 

/3o2 

X 

Z 

/3oi 

/3o2 

X 

Z 

/3oi 

/3o2 

X 


ind 

1.0 

-0.2 

-1.0 

0.3 

1.00 

1.00 

Complete 

1.00 1.00 

0.96 

0.96 

0.95 

0.93 

exch 

1.2 

-0.1 

0.1 

0.2 

1.00 

1.00 

1.00 

1.26 

0.97 

0.96 

0.96 

0.95 

unst 

0.2 

-0.2 

-1.0 

0.3 

1.00 

1.00 

1.00 

1.26 

0.96 

0.95 

0.96 

0.95 

unit 

0.9 

-0.2 

-0.9 

0.2 

1.00 

1.01 

1.00 

1.28 

0.96 

0.95 

0.95 

0.94 

cat. exch 

0.9 

-0.2 

-1.0 

-0.2 

1.00 

1.01 

1.00 

1.28 

0.96 

0.95 

0.95 

0.94 

time.exch 

0.4 

0.0 

0.5 

0.2 

1.00 

1.00 

1.00 

1.28 

0.96 

0.95 

0.94 

0.95 

RC 

0.0 

0.4 

1.4 

0.1 

1.00 

1.00 

1.00 

1.28 

0.95 

0.95 

0.96 

0.96 

ind 

2.2 

-0.6 

-0.3 

0.6 

1.00 

1.00 

WGEE(r+) 

1.00 1.00 

0.98 

0.97 

0.97 

0.95 

exch 

1.7 

-0.3 

0.9 

-0.3 

0.99 

0.98 

0.99 

1.20 

0.97 

0.96 

0.95 

0.96 

unst 

2.2 

-0.8 

-1.5 

-0.8 

0.99 

0.98 

0.99 

1.20 

0.97 

0.96 

0.97 

0.97 

unit 

2.0 

-0.8 

-0.7 

-1.7 

0.99 

0.99 

0.99 

1.20 

0.97 

0.96 

0.96 

0.95 

cat. exch 

2.1 

-0.9 

-0.7 

-2.1 

0.99 

0.99 

0.99 

1.20 

0.97 

0.96 

0.96 

0.95 

time.exch 

2.3 

-0.7 

-0.9 

-1.5 

0.99 

0.99 

0.99 

1.20 

0.96 

0.95 

0.95 

0.95 

RC 

1.7 

-0.1 

1.5 

0.3 

0.99 

0.99 

0.99 

1.20 

0.97 

0.96 

0.96 

0.95 

ind 

3.6 

-1.1 

-5.6 

-0.2 

1.00 

MIGEE(a:+) 

1.00 1.00 1.00 

0.96 

0.95 

0.97 

0.92 

exch 

4.3 

-1.1 

-5.0 

-0.3 

0.99 

0.99 

1.00 

1.21 

0.95 

0.95 

0.96 

0.93 

unst 

4.2 

-1.5 

-7.5 

-0.9 

0.99 

0.99 

1.00 

1.21 

0.95 

0.95 

0.97 

0.95 

unit 

3.6 

-1.2 

-5.6 

-0.8 

1.00 

1.00 

1.01 

1.22 

0.96 

0.96 

0.96 

0.94 

cat. exch 

3.7 

-1.3 

-5.8 

-1.4 

1.00 

1.00 

1.01 

1.22 

0.96 

0.96 

0.96 

0.94 

time.exch 

3.9 

-1.4 

-6.0 

-1.0 

1.00 

1.00 

1.01 

1.22 

0.95 

0.94 

0.95 

0.94 

RC 

3.4 

-0.8 

-4.0 

-0.8 

1.00 

1.00 

1.01 

1.22 

0.94 

0.96 

0.96 

0.95 

ind 

-0.3 

0.2 

-0.5 

1.0 

1.00 

DRGEE(x+ 
1.00 1.00 

,r+) 

1.00 

0.97 

0.96 

0.97 

0.94 

exch 

-0.6 

0.3 

0.2 

0.9 

1.01 

0.98 

0.99 

1.17 

0.96 

0.94 

0.96 

0.93 

unst 

-0.3 

-0.1 

-1.8 

0.5 

1.01 

0.98 

0.99 

1.16 

0.95 

0.94 

0.96 

0.95 

unit 

-0.5 

0.0 

-0.9 

0.1 

1.00 

1.00 

1.00 

1.17 

0.95 

0.95 

0.96 

0.95 

cat. exch 

-0.5 

0.0 

-1.0 

-0.3 

1.00 

1.00 

1.00 

1.17 

0.95 

0.95 

0.96 

0.95 

time.exch 

0.4 

0.0 

-1.6 

0.4 

1.00 

1.00 

0.99 

1.17 

0.94 

0.95 

0.95 

0.95 

RC 

-1.0 

0.7 

1.4 

1.6 

1.00 

1.00 

1.00 

1.16 

0.95 

0.95 

0.95 

0.95 

ind 

0.9 

-0.3 

-3.6 

0.6 

1.00 

DRGEE(a;- 
1.00 1.00 

, r+) 
1.00 

0.97 

0.96 

0.97 

0.93 

exch 

0.4 

0.0 

-2.3 

0.9 

0.99 

0.98 

1.00 

1.18 

0.96 

0.95 

0.95 

0.92 

unst 

0.9 

-0.5 

-4.7 

0.4 

0.99 

0.98 

1.01 

1.18 

0.96 

0.96 

0.96 

0.95 

unit 

0.7 

-0.4 

-4.0 

0.2 

1.00 

1.00 

1.01 

1.18 

0.96 

0.96 

0.96 

0.95 

cat. exch 

0.8 

-0.4 

-4.1 

-0.2 

1.00 

1.00 

1.01 

1.19 

0.96 

0.96 

0.96 

0.95 

time.exch 

1.2 

-0.2 

-3.7 

0.7 

0.99 

0.99 

1.00 

1.18 

0.95 

0.95 

0.96 

0.95 

RC 

0.3 

0.3 

-1.9 

1.8 

1.00 

0.99 

1.01 

1.18 

0.96 

0.96 

0.96 

0.95 

ind 

0.0 

0.0 

-1.4 

0.7 

1.00 

DRGEE(a;+ 
1.00 1.00 

,r-) 

1.00 

0.97 

0.96 

0.96 

0.94 

exch 

-0.2 

0.2 

-0.6 

0.1 

0.99 

0.98 

0.99 

1.15 

0.95 

0.94 

0.95 

0.93 

unst 

-0.1 

-0.1 

-2.1 

-1.0 

0.99 

0.98 

0.99 

1.15 

0.95 

0.94 

0.96 

0.95 

unit 

-0.2 

-0.1 

-1.9 

-0.3 

1.00 

0.99 

0.99 

1.16 

0.95 

0.95 

0.96 

0.95 

cat. exch 

-0.1 

-0.2 

-1.9 

-0.7 

1.00 

0.99 

0.99 

1.16 

0.95 

0.95 

0.96 

0.95 

time.exch 

0.6 

-0.2 

-2.8 

-0.5 

0.99 

0.99 

0.99 

1.16 

0.94 

0.94 

0.94 

0.95 

RC 

-0.5 

0.5 

0.1 

1.1 

1.00 

0.99 

0.99 

1.15 

0.94 

0.94 

0.96 

0.95 


indicates correctly specified model and “ ” indicates misspecified model omitting the Z\ predictor 





Table 10: Convergence Rate (CR) for the Simulation Study 




sample size 


50 

150 

300 

600 

ind 

0.97 

1.00 

1.00 

1.00 

exch 

0.59 

0.89 

0.96 

1.00 

unst 

0.08 

0.59 

0.89 

0.99 

unif 

0.97 

1.00 

1.00 

1.00 

cat. exch 

0.97 

1.00 

1.00 

1.00 

time.exch 

0.96 

1.00 

1.00 

1.00 

RC 

0.67 

0.99 

1.00 

1.00 


Table [10] shows the convergence rate obtained from the simulation results for seven 
association structures and four sample sizes. The convergence rate for the working asso¬ 
ciation structure (C) was calculated as 


CR^ 


1000 

total number of simulations 


The local odds ratio parametrization presents a clear advantage over the correlation coef- 
hcient in terms of convergence issues, specially for unstructured association matrices and 
small sample sizes. 
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A.3 Additional Real Data Results 
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Table 11: Results for Rheumatic Mitral Stenosis study under exchangeability, time ex¬ 
changeability and RC association structures 



Exchangeable 

Time Exchangeability 


RC 


Parameter 

Est. 

SE 

P 

Est. 

SE 

P 

Est. 

SE 

P 


-1.064 

0.782 

0.174 

-0.909 

Available 

0.730 

0.213 

-0.881 

0.744 

0.236 

1^02 

0.721 

0.776 

0.353 

0.890 

0.724 

0.219 

0.924 

0.738 

0.210 

Success 

0.680 

0.397 

0.086 

0.626 

0.387 

0.105 

0.641 

0.391 

0.101 

Total Score 

-0.108 

0.107 

0.313 

-0.130 

0.101 

0.197 

-0.137 

0.102 

0.179 

Cn 

0.468 

0.339 

0.168 

0.477 

0.314 

0.129 

0.496 

0.317 

0.118 

time=2 

1.083 

0.398 

0.007 

1.022 

0.383 

0.008 

1.081 

0.390 

0.006 

time=3 

1.208 

0.402 

0.003 

1.120 

0.380 

0.003 

1.214 

0.392 

0.002 

Poi 

-1178 

0.807 

0.144 

-0.934 

WGEE 

0.756 

0.217 

-0.892 

0.771 

0.248 

1^02 

0.611 

0.797 

0.443 

0.848 

0.741 

0.252 

0.898 

0.757 

0.236 

Success 

0.744 

0.410 

0.070 

0.623 

0.400 

0.119 

0.645 

0.404 

0.111 

Total Score 

-0.101 

0.109 

0.355 

-0.127 

0.103 

0.216 

-0.136 

0.105 

0.195 

Cn 

0.490 

0.351 

0.162 

0.433 

0.322 

0.179 

0.449 

0.325 

0.167 

time=2 

1.037 

0.414 

0.012 

1.013 

0.398 

0.011 

1.060 

0.404 

0.009 

time=3 

1.395 

0.418 

0.001 

1.363 

0.405 

0.001 

1.437 

0.415 

0.001 

Poi 

-1.100 

0.618 

0.075 

-0.865 

MIGEE 

0.573 

0.131 

-0.908 

0.567 

0.109 

(^02 

0.743 

0.621 

0.231 

0.973 

0.575 

0.091 

0.934 

0.572 

0.102 

Success 

1.079 

0.334 

0.001 

0.834 

0.322 

0.010 

0.887 

0.324 

0.006 

Total Score 

-0.098 

0.083 

0.239 

-0.131 

0.078 

0.094 

-0.131 

0.078 

0.093 

Cn 

0.237 

0.280 

0.397 

0.281 

0.265 

0.288 

0.350 

0.272 

0.198 

time=2 

1.095 

0.317 

0.001 

1.075 

0.320 

0.001 

1.033 

0.323 

0.001 

time=3 

1.071 

0.316 

0.001 

1.091 

0.308 

0.000 

1.072 

0.313 

0.001 

f^Ol 

-1.199 

0.781 

0.125 

-0.896 

DRGEE 

0.692 

0.195 

-0.852 

0.700 

0.224 

1^02 

0.591 

0.785 

0.451 

0.888 

0.686 

0.195 

0.934 

0.696 

0.180 

Success 

0.808 

0.380 

0.034 

0.727 

0.368 

0.048 

0.748 

0.371 

0.044 

Total Score 

-0.097 

0.108 

0.368 

-0.137 

0.094 

0.148 

-0.143 

0.096 

0.135 

Cn 

0.475 

0.373 

0.202 

0.490 

0.326 

0.133 

0.495 

0.330 

0.134 

time=2 

0.937 

0.630 

0.137 

0.985 

0.534 

0.065 

0.983 

0.540 

0.069 

time=3 

1.203 

0.601 

0.046 

1.189 

0.498 

0.017 

1.208 

0.507 

0.017 


31 






